Civil Engineering Reference
In-Depth Information
f1=(c1+pt5*c2*fm)*mv+(one+pt5*c2*fk)*kv; CALL sparin(f1,kdiag); time=zero
!-----------------------time stepping loop--------------------------------
WRITE(11,'(/A,I5))')" Result at node",nres
WRITE(11,'(A)')"
time
load
x-disp
y-disp"
WRITE(11,'(4E12.4)')time,load(time),x0(nf(:,nres))
timesteps: DO j=1,nstep
time=time+dtim; loads=zero
x1=(c1+pt5*c2*fm)*x0+(c2+two*fm)*d1x0+(two+c4*fm)*d2x0
DO i=1,loaded_nodes
loads(nf(:,node(i)))=
&
val(i,:)*(theta*load(time)+(one-theta)*load(time-dtim))
END DO; CALL linmul_sky(mv,x1,d1x1,kdiag); d1x1=loads+d1x1
loads=pt5*c2*fk*x0+two*fk*d1x0+c4*fk*d2x0
CALL linmul_sky(kv,loads,x1,kdiag); x1=x1+d1x1; CALL spabac(f1,x1,kdiag)
d2x1=d2x0+((x1-x0)*c1-d1x0*c2-d2x0*two-d2x0)/theta
d1x1=d1x0+pt5*dtim*(d2x1+d2x0); x0=x0+dtim*d1x0+two*c3*d2x0+d2x1*c3
d1x0=d1x1; d2x0=d2x1
IF(j/npri*npri==j)WRITE(11,'(4E12.4)')time,load(time),x0(nf(:,nres))
END DO timesteps
STOP
CONTAINS
FUNCTION load(t) RESULT(load_result)
!-----------------------Load-time function--------------------------------
IMPLICIT NONE
INTEGER,PARAMETER::iwp=SELECTED_REAL_KIND(15)
REAL(iwp),INTENT(IN)::t
REAL(iwp)::load_result
load_result=COS(0.3_iwp*t)
RETURN
END FUNCTION load
END PROGRAM p114
New scalar reals:
d6
set to 6.0
two
set to 2.0
This algorithm is described in Section 3.13.3. The essential step is that shown in (3.142),
which is of exactly the same form as (3.139) for the “theta” method, and so this program can
be expected to resemble the previous one very closely. The structure chart of Figure 11.8
is again appropriate and the problem layout and data given in Figure 11.11 are essentially
the same as in Figure 11.9. No new variables are involved, but the parameter
θ
is set to
4.0
16
1
4
1.0
2
17
18
5
3
E= 1 kN/m
2
n
= 0.3
r
= 1 t/m
3
F= cos
w
t
Figure 11.11
Mesh and data for Program 11.4 example (
Continued on page 486
)
