[分享]精细时程积分FORTRAN源程序

风花雪月 · 发表于 2005-8-7 16:20

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。

您需要登录才可以下载或查看，没有账号？我要加入

x

主程序main.f90，搞了两个例子进行调用，例子是钟院士那篇文献里的，附后的文献里有的。一个程序是reciseTIM.f90，里面包含了精细时程积分的子程序。

最后的PDF文件是程序说明.

主程序Program Main

use PreciseTimeIntegration
Real(8), Dimension(8, 8) :: K,M,C,INVK
Real(8), Dimension(8, 4) :: R,U,V,A
Real(8) dt
dt=1
R=0
C=0
C(7,7)=5
C(7,8)=-5
C(8,7)=-5
C(8,8)=5
DO i=1,8
M(i,i)=8
END DO
DO i=1,7
K(i,i+1)=-4
k(i+1,i)=-4
END DO
DO i=2,7
k(i,i)=8
END DO
K(1,1)=4
K(8,8)=4
! WRITE (*,*) M
! WRITE (*,*) K
! WRITE (*,*) C
! INVK=INV(K,8)
! STOP
U=0
U(8,1)=10
V=0
A=0
call PreciseTIM (M, C, k, R, 8, 4, dt, U, V, A)
write (*,*) U
! write (*,*) V
! write (*,*) A
end

复制代码

本程序来自http://bbs.tongji.edu.cn/bbsanc.php?path=%2FPersonalCorpus%2FA%2FABAYA%2FD4F47E4DC%2FA41E7F755

风花雪月 · 发表于 2005-8-7 16:21

精细时程积分文件：

PreciseTIM.f90

!-----------------------------------------------------------------
Module PreciseTimeIntegration

Contains

!//////////////////////////////////////////////////////////////
Subroutine PreciseTIM (M, G, K, R, N, Nt, tao, X, XX, XXX)
!//////////////////////////////////////////////////////////////
! 精细时程积分法程序，
! ABAYA，同济大学建筑工程系，2001.7.20
!!//////////////////////////////////////////////////////////////
Real*8 M(N,N), G(N,N), K(N,N), R(N,Nt), X(N,Nt), &
XX(N,Nt), XXX(N,Nt), H(2*N,2*N), T(2*N,2*N), &
R0(2*N,Nt), R1(2*N,Nt), B(N,N), C(N,N), tao
Call CalH (M, G, K, N, Nt, H, B, C)
Call CalT (H, tao, N, T)
Call CalR0R1 (R, R0, R1, N, Nt)
Call CalX (T, H, R, R0, R1, tao, N, Nt, M, G, B, C, X, XX, XXX)
End Subroutine PreciseTIM
!//////////////////////////////////////////////////////////////
Subroutine CalH (M, G, K, N, Nt, H, B, C)
Real*8 M(N,N), G(N,N), K(N,N), R(N,Nt), INVM(N,N), &
A(N,N), B(N,N), C(N,N), D(N,N), H(2*N,2*N)
INVM=INV(M,N)
A=-0.5*Matmul(INVM,G)
B=0.25*Matmul(Matmul(G,INVM),G)-K
C=-0.5*Matmul(G,INVM)
D=INVM
H(1:N,1:N)=A
H(N+1:2*N,1:N)=B
H(1:N,N+1:2*N)=D
H(N+1:2*N,N+1:2*N)=C
End Subroutine CalH
!//////////////////////////////////////////////////////////////
Subroutine CalT (H, tao, N, T)
Real*8, Dimension(2*N,2*N) :: H, T, Ta, I
Real*8 tao, dt
Integer m
m=2**20
dt=tao/m
I=0; T=0
Do j=1,2*N
I(j,j)=1
End do
Ta=Matmul(dt*H,(I+0.5*dt*H))
Do j=1,20
Ta=2*Ta+Matmul(Ta,Ta)
End do
T=I+Ta
End Subroutine CalT
!//////////////////////////////////////////////////////////////
Subroutine CalR0R1 (R, R0, R1, N, Nt)
Real*8 R(N,Nt), R0(2*N,Nt), R1(2*N,Nt)
R0=0; R1=0;
R0(N+1:2*N,:)=R
Do i=2,Nt
R1(N+1:2*N,i-1)=R(:,i)-R(:,i-1)
End do
End Subroutine CalR0R1
!//////////////////////////////////////////////////////////////
Subroutine CalX (T, H, R, R0, R1, tao, N, Nt, M, G, B, C, X, XX, XXX)
Real*8 X(N,Nt), XX(N,Nt), XXX(N,Nt), H(2*N,2*N), &
T(2*N,2*N), R0(2*N,Nt), R1(2*N,Nt), V(2*N,Nt), &
p(N,Nt), q(N,Nt), R(N,Nt), B(N,N), C(N,N), M(N,N), &
G(N,N), tao, INVH(2*N,2*N), index
q=X
p=Matmul(M,XX)+0.5*Matmul(G,X)
V(1:N,:)=q
V(N+1:2*N,:)=p
index=0
Do i=1,2*N
Do j=1,Nt
IF (ABS(R0(i,j)).GT.1E-8) THEN
index=1; GOTO 10
END IF
End do
End do
10 If (abs(index).gt.1e-8) then
INVH=INV(H,2*N)
Do i=2,Nt
WRITE (*,'("*********** LOAD STEP: ",I5," ***********")') i
V(:,i)=Matmul(T,(V(:,i-1)+Matmul(INVH,(R0(:,i-1)+ &
Matmul(INVH,R1(:,i-1))))))-Matmul(INVH &
,(R0(:,i-1)+Matmul(INVH,R1(:,i-1))+ &
tao*R1(:,i-1)))
End do
Else
Do i=2,Nt
WRITE (*,'("*********** LOAD STEP: ",I5," ***********")') i
V(:,i)=Matmul(T,V(:,i-1))
End do
End if
q=V(1:N,:)
p=V(N+1:2*N,:)
Do i=1,Nt
X(:,i)=q(:,i)
XX(:,i)=Matmul(INV(M,N),p(:,i))-0.5*Matmul(Matmul(INV(M,N),G),X(:,i))
XXX(:,i)=Matmul(Matmul(INV(M,N),B),X(:,i))-0.5*Matmul(Matmul(INV(M,N) &
,G),XX(:,i))+Matmul(Matmul(INV(M,N),C),p(:,i))+Matmul(INV(M,N),R(:,i))
End do
End Subroutine CalX
!//////////////////////////////////////////////////////////////
Function INV (A, N)
Use Numerical_Libraries
REAL*8, DIMENSION(N,N) :: A, INV
INV=0
CALL DLINRG (N, A, N, INV, N)
End Function INV
!//////////////////////////////////////////////////////////////

复制代码

End Module PreciseTimeIntegration

linqus · 发表于 2005-8-7 18:54

赞一个先,
自己也曾编过线性的wilson时程分析程序.
可以解决剪切层模型问题.没有作验证:)

wgp1976 · 发表于 2005-8-8 15:35

对此比较感兴趣，先睹为快

wangjiaqi4 · 发表于 2005-8-9 14:27

谢谢了，感谢共享这么好的源代码。

simon21 · 发表于 2005-8-9 21:23

谢谢共享

mxtfirst · 发表于 2005-8-16 10:37

pdf的文档怎么无法下载？

bgzhoumse · 发表于 2005-8-22 20:39

使用marc进行有限元分析时一定要用fortran语言进行有关的编程吗？可不可以用C语言

NASA · 发表于 2005-8-22 20:51

Marc的接口是Fortran的，会C语言学Fortran不难
Marc各个版本与fortran版本搭配如下：

marc2003需要fortran6.0以上版本
marc2000以及marc2001需要fortran5.0以上版本

dovenkai · 发表于 2005-8-23 20:26

ding

liu3806881 · 发表于 2005-8-27 22:39

哈哈，找拉好久的拉
好！！！！

xj2070 · 发表于 2005-8-28 16:09

非常不错的东西，感谢搂住的热心

hit_zzc · 发表于 2005-8-30 17:22

好咚咚

johhan · 发表于 2005-8-30 20:40

先看一下,是否有用

jsswdd · 发表于 2005-9-17 01:34

xiexie

账号		自动登录	找回密码
密码			我要加入

[Fortran] [分享]精细时程积分FORTRAN源程序

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。