<边界元法计算界面几何性质 >的数据输入问题
!***************************************************************!
! CALCULATION OF THE GEOMETRICAL PROPERTIES OF SECTION (GPS)
! (boundary integration methon)
! SELINA999 2023-06-04
!***************************************************************
!
IMPLICIT REAL*8(A-H,O-Z)
CHARACTER*8 DATA
CHARACTER TITLE(70)
DIMENSION GPS(6),RES(6),ROC(101),X(101),Y(101)
!
! NOA------------total number of part areas
! SOA------------sign of the part area
! NP-------------total number of nodes of the part area
! X(101),Y(101)--node coordinates
! NE-------------total number of boundary elements of the
! area
! ROC(101)-------radiu of the boundary element (POSITIVE
! radiu for convex boundary, NEGATIVE radiu
! for concave boundary,'0'for straight line)
! NDV------------number of divisions of each boundary element
! integration
! RES(6),GPS(6)--GPS relative to global axis of the part and
! total area respectively,(1)A,(2)Sx,(3)Sy,
! (4)Ix,(5)Iy,(6)Ixy
! XC,YC----------centriod coordinates
! AIX,AIY,AIXY---Ix,Iy,Ixy relative to parallel centriod axis
! PIMAX,PIMIN----centriodal principle moments of inertia
! ANG------------declination of the principle axis relative to
! Imax
!
WRITE(*,*) 'Input the file name of data please'
READ(*,'(A)') DATA
OPEN(5,FILE=DATA)
OPEN(6,FILE='RESU',STATUS='NEW')
!
! a.Input the governing data
!
READ(5,2000) TITLE
WRITE(6,2010) TITLE
WRITE(*,2010) TITLE
!
WRITE(*,2020)
READ(*,*) NDV
READ(5,*) NOA
WRITE(6,2030) NOA
WRITE(*,2030) NOA
!
DO 10 I=1,6
10 GPS(I)=0.0
!
DO 90 M=1,NOA
!
! Input part area data
!
READ(5,*) SOA
READ(5,*) NP
WRITE(6,2040) M,SOA,NP
WRITE(*,2040) M,SOA,NP
DO 15 I=1,NP
READ(5,*) X(I),Y(I),ROC(I)
WRITE(6,2050) I,X(I),Y(I),I,ROC(I)
15 WRITE(*,2050) I,X(I),Y(I),I,ROC(I)
X(NP+1)=X(1)
Y(NP+1)=Y(1)
!
DO 20 I=1,6
20 RES(I)=0.0
!
NE=NP
DO 70 N=1,NE
!
! b.Determime the element constants
!
DX=X(N+1)-X(N)
DY=Y(N+1)-Y(N)
CL=DSQRT(DX*DX+DY*DY)
!
IF (ROC(N).EQ.0.) GOTO 30
!
! For curve boundary element
!
DCX=DY/CL
DCY=-DX/CL
XM=0.5*(X(N+1)+X(N))
YM=0.5*(Y(N+1)+Y(N))
R=DABS(ROC(N))
S=ROC(N)/R
D2=R*R-CL*CL/4.
D=DSQRT(D2)
XCC=XM-S*D*DCX
YCC=YM-S*D*DCY
ZETA=2.0*DASIN(0.5*CL/R)
ALFA=DACOS(DCX)
IF (DCY.GE.0.) GOTO 30
ALFA=-ALFA
!
30 X1=X(N)
Y1=Y(N)
!
! c.Integrate along the boundary element
!
DO 60 I=1,NDV
PDV=1.0/NDV
!
IF (ROC(N).NE.0.) GOTO 40
X2=X1+PDV*DX
Y2=Y1+PDV*DY
GOTO 50
!
40 BI=I
BETA=ALFA-S*0.5*ZETA+S*BI*PDV*ZETA
X2=XCC+S*R*DCOS(BETA)
Y2=YCC+S*R*DSIN(BETA)
!
! Detetmine the integral point,length and direction cosine
!
50 DDX=X2-X1
DDY=Y2-Y1
DL=DSQRT(DDX*DDX+DDY*DDY)
DDCX=DDY/DL
DDCY=-DDX/DL
!
DXM=0.5*(X1+X2)
DYM=0.5*(Y1+Y2)
RES(1)=RES(1)+DXM*DDCX*DL
RES(2)=RES(2)+0.5*DYM*DYM*DDCY*DL
RES(3)=RES(3)+0.5*DXM*DXM*DDCX*DL
RES(4)=RES(4)+DYM**3*DDCY*DL/3.
RES(5)=RES(5)+DXM**3*DDCX*DL/3.
RES(6)=RES(6)+0.5*DXM*DXM*DYM*DDCX*DL
X1=X2
Y1=Y2
60 CONTINUE
70 CONTINUE
!
! Detemine GPS of total area
!
DO 80 I=1,6
80 GPS(I)=GPS(I)+RES(I)*SOA
90 CONTINUE
!
WRITE(6,2060) NDV
WRITE(*,2060) NDV
!
! d.Calculate GPS respect to parallel central axis
!
A=GPS(1)
XC=GPS(3)/A
YC=GPS(2)/A
AIX=GPS(4)-YC*YC*A
AIY=GPS(5)-XC*XC*A
AIXY=GPS(6)-XC*YC*A
!
! e.Calculate the centriodal principle moments of inertia
!
C1=0.5*(AIX+AIY)
C2=0.5*(AIX-AIY)
C3=DSQRT(C2*C2+AIXY*AIXY)
PIMAX=C1+C3
PIMIN=C1-C3
IF (AIX.EQ.PIMIN) GOTO 95
ANG=DATAN(-AIXY/(AIX-PIMIN))*180.0/3.1415926
GOTO 100
95 ANG=-90.0
100 WRITE(6,2070) A,XC,YC,AIX,AIY,AIXY,PIMAX,PIMIN,ANG
WRITE(*,2070) A,XC,YC,AIX,AIY,AIXY,PIMAX,PIMIN,ANG
!
!
2000 FORMAT(70A1)
2010 FORMAT(2X,70A1)
2020 FORMAT(/,2X,'please input the total number of divisions of &
& each bounday element for integrationNDV')
2030 FORMAT(/,2X,'TOTAL NUMBER OF PART AREAS',6X,'NOA=',I4)
2040 FORMAT(//2X,'PART AREA NO:',I3,/,8X,'SIGN OF THE PART AREA &
& SOA=',2X,F3.0,/,8X,'TOTAL NUMBER NODES NP=',I4,//,&
& 2X,'NODE',8X,'COORDINATES(MM)',8X,'ELEMENT',5X,'RADIU(MM)',/,&
& 4X,'NO',10X,'X',9X,'Y',13X,'NO',10X,'R',/)
2050 FORMAT(3X,I3,5X,2F10.2,8X,I3,5X,F10.2)
2060 FORMAT(//,2X,'GEOMETRICAL PROPERTIES OF THE SECTION (NDV=',&
& I4,')',//)
2070 FORMAT(4X,'TOTAL AREA',22X,'A=',E11.4,2X,'(MM**2)',//,&
& 4X,'CENTRION COORDINATES',11X,'Xc=',E11.4,2X,'(MM)',/,&
& 35X,'Yc=',E11.4,2X,'(MM)',//&
& 4X,'MOMENT OF INTRIA RELATIVE TO THE PARALELL CENTRIOD AXIS',/,&
&/,35X,'Ix=',E11.4,2X,'(MM**4)',/,35X,'Iy=',E11.4,2X,'(MM**4)',&
&/,34X,'Ixy=',E11.4,2X,'(MM**4)',//,&
& 4X,'CENTRIODAL PRINCIPLE MOMENTS OF INERTIA',/,&
&/,33X,'Imax=',E11.4,2X,'(MM**4)',/,33X,'Imin=',E11.4,2X,'(MM**4)' &
&,//,4X,'DECLINATION OF THE PRINCIPLE AXIS OF Imax',/,&
&/,32X,'ANGLE=',E11.4,2X,'(DEG)')
STOP
END
书上的例题：
RECTANGLE WITH TWO SEMI-CICLE CUT
1
1.
8
0.,0.,0.
80.,0.,0.
80.,18.,-32.
80.,82.,0.
80.,100.,0.
0.,100.,0.
0.,82.,-32.
0.,18.,0.
编译可以通过，但数据好像哪里输入不对？
我在提示输入NDV时输入1，按你的输入数据建立了一个文件，程序执行正常结束 necrohan 发表于 2023-6-12 12:32
我在提示输入NDV时输入1，按你的输入数据建立了一个文件，程序执行正常结束 ...
能不能发个编译好的程序给我，谢谢！ 我编译完成，输入文件名后，直接就退出来了，真么回事？:-sleepy: test_read
https://www.aliyundrive.com/s/qaL2oyguTHy
提取码: 53be
点击链接保存，或者复制本段内容，打开「阿里云盘」APP ，无需下载极速在线查看，视频原画倍速播放。 necrohan 发表于 2023-6-13 09:13
test_read
https://www.aliyundrive.com/s/qaL2oyguTHy
提取码: 53be
我编译完成，运行的时候出错:-o c:went3 本帖最后由 selina999 于 2023-6-13 14:46 编辑
necrohan 发表于 2023-6-13 09:13
test_read
https://www.aliyundrive.com/s/qaL2oyguTHy
提取码: 53be
谢谢，按照您发过来的可以运行，但是我编译的就出错，哪的问题:-loveliness: 本帖最后由 necrohan 于 2023-6-14 08:34 编辑
看错误提示，GPS2.f90文件的第37行出错，因为文件resu.txt已经存在所以无法打开。
不知道你要干什么，如果这个文件只是输出用可以在打开文件的语句直接覆盖，一般文件句柄不要用5、6，那是系统显示用的。
先把resu.txt删除再运行试试。
页:
[1]