Civil Engineering Reference
In-Depth Information
Program 12.9 Forced vibration analysis of a three-dimensional elastic solid. Implicit
integration in time. Compare Program 11.4.
PROGRAM p129
!-------------------------------------------------------------------------
! Program 11.6 forced vibration of a 3 - d elastic
! solid using uniform 20-node hexahedral elements (nxe even)
! numbered in the x-z direction - lumped or consistent mass
! implicit integration by theta method : parallel version
!-------------------------------------------------------------------------
USE new_library; USE geometry_lib ; USE precision; USE gather_scatter6
USE global_variables1;USE timing;USE mp_module;USE utility; IMPLICIT NONE
INTEGER::nxe,nye,nze,nn,nr,nip,nodof=3,nod=20,nst=6,neq_temp,nn_temp,
&
i,j,k,iel,ndim=3,nstep,npri,iters,limit , it,is ,nres
! ndof, nels, ntot, neq are now global variables ; not declared
REAL(iwp)::aa,bb,cc,e,v,det,rho,alpha1,beta1,omega,theta,period,pi,dtim, &
volume,c1,c2,c3,c4,real_time,tol,big,up,alpha,beta
CHARACTER(LEN=15)::element='hexahedron'
LOGICAL :: consistent = .FALSE. , converged
!----------------------------- dynamic arrays-----------------------------
REAL(iwp),ALLOCATABLE::loads_pp(:),points(:,:),dee(:,:),coord(:,:), &
fun(:),jac(:,:), der(:,:),deriv(:,:), weights(:), &
bee(:,:),km(:,:),p_g_co_pp(:,:,:),x1_pp(:), &
d1x1_pp(:),d2x1_pp(:),emm(:,:),ecm(:,:),x0_pp(:), &
d1x0_pp(:),d2x0_pp(:),store_km_pp(:,:,:),vu_pp(:), &
store_mm_pp(:,:,:),u_pp(:),p_pp(:),d_pp(:), &
x_pp(:),xnew_pp(:),pmul_pp(:,:),utemp_pp(:,:), &
diag_precon_pp(:),diag_precon_tmp(:,:),temp_pp(:,:,:)
INTEGER, ALLOCATABLE::rest(:,:), g(:), num(:), g_num_pp(:,:), g_g_pp(:,:)
!------------------------input and initialisation-------------------------
timest(1) = elap_time( ) ;
CALL find_pe_procs(numpe,npes)
IF(numpe==npes) THEN
OPEN (10,FILE='p129.dat',STATUS= 'OLD',ACTION='READ')
READ (10,*) nels,nxe,nze,nip,aa,bb,cc,rho,e,v,
&
alpha1,beta1,nstep,npri,theta,omega,tol,limit
END IF
CALL bcast_inputdata_p129(numpe,npes,nels,nxe,nze,nip,aa,bb,cc,rho,e,v, &
alpha1,beta1,nstep,theta,npri,omega,tol,limit)
CALL calc_nels_pp ; ndof=nod*nodof ; neq_temp = 0 ; nn_temp = 0
ntot = ndof; nye = nels/nxe/nze; nr=3*nxe*nze+2*nxe+2*nze+1
nres = 3*(nye*(nxe+1)*(nze+1)+nr*(nye-1)+(nxe+1))
ALLOCATE(rest(nr,nodof+1),points(nip,ndim),g(ntot),fun(nod), &
dee(nst,nst),coord(nod,ndim),jac(ndim,ndim),weights(nip), &
der(ndim,nod), deriv(ndim,nod), bee(nst,ntot), km(ntot,ntot), &
num(nod),g_num_pp(nod,nels_pp),g_g_pp(ntot,nels_pp),
&
emm(ntot,ntot),ecm(ntot,ntot),p_g_co_pp(nod,ndim,nels_pp),
&
store_km_pp(ntot,ntot,nels_pp),utemp_pp(ntot,nels_pp),
&
pmul_pp(ntot,nels_pp),store_mm_pp(ntot,ntot,nels_pp),
&
temp_pp(ntot,ntot,nels_pp),diag_precon_tmp(ntot,nels_pp))
rest = 0; DO i=1,nr; rest(i,1) = i; END DO ; ielpe = iel_start
pi=ACOS(-1._iwp) ; period = 2._iwp*pi/omega ; dtim =period/20._iwp
c1=(1._iwp-theta)*dtim; c2=beta1-c1
c3=alpha1+1._iwp/(theta*dtim); c4=beta1+theta*dtim
CALL deemat (dee,e,v); CALL sample(element,points,weights)
!-------------- loop the elements to set up global arrays ----------------
elements_1: DO iel = 1 , nels_pp
CALL geometry_20bxz(ielpe,nxe,nze,aa,bb,cc,coord,num)
CALL find_g( num , g , rest ) ; g_num_pp(:,iel) = num
