Civil Engineering Reference
In-Depth Information
!-------------------- get starting r--------------------------------------
IF(loaded_freedoms>0) THEN
CALL reindex_fixed_nodes
&
(ieq_start,no,no_local_temp,num_no,no_index_start)
ALLOCATE(no_local(1:num_no)) ; no_local = no_local_temp(1:num_no)
DEALLOCATE(no_local_temp)
DO i = 1 , num_no
r_pp(no_local(i)-ieq_start+1) = val(no_index_start + i - 1)
END DO
END IF ; q = SUM_P(r_pp)
IF(numpe==1) WRITE(11,'(A,E12.4)') "The total load is ", q
diag_precon_pp=1._iwp/diag_precon_pp;d_pp=diag_precon_pp*r_pp;p_pp=d_pp
!--------------------preconditioned c. g. iterations----------------------
iters = 0
iterations : DO
iters=iters+1; u_pp=0._iwp;pmul_pp=.0_iwp;CALL gather(p_pp,pmul_pp)
elements_2 : DO iel = 1, nels_pp
utemp_pp(:,iel) = MATMUL(km,pmul_pp(:,iel))
!
CALL dgemv('n',ntot,ntot,1.0,km,ntot,
&
! pmul_pp(:,iel),1,0.0,utemp_pp(:,iel),1)
END DO elements_2 ; CALL scatter(u_pp,utemp_pp)
!--------------------------pcg equation solution--------------------------
up=DOT_PRODUCT_P(r_pp,d_pp);alpha= up/ DOT_PRODUCT_P(p_pp,u_pp)
xnew_pp = x_pp + p_pp* alpha ; r_pp=r_pp - u_pp*alpha
d_pp = diag_precon_pp*r_pp; beta=DOT_PRODUCT_P(r_pp,d_pp)/up
p_pp=d_pp+p_pp*beta
CALL checon_par(xnew_pp,x_pp,tol,converged,neq_pp)
IF(converged .OR. iters==limit) EXIT
END DO iterations
IF(numpe==1)THEN
WRITE(11,'(A,I5)')"The number of iterations to convergence was ",iters
WRITE(11,'(A,E12.4)')"The central nodal displacement is :",xnew_pp(1)
END IF
!-------------------recover stresses at centroidal gauss-point------------
nip=1; points = .0_iwp ; eld_pp = .0_iwp; CALL gather(xnew_pp,eld_pp)
iel = 1; coord= p_g_co_pp(:,:,iel); eld=eld_pp(:,iel)
IF(numpe==1)WRITE(11,'(A)')"The Centroid point stresses for element 1 are"
gauss_pts_2: DO i= 1 , nip
CALL shape_der(der,points,i); jac= MATMUL(der,coord)
CALL invert (jac); deriv= MATMUL(jac,der) ; CALL beemat(bee,deriv)
eps = MATMUL (bee,eld) ; sigma = MATMUL (dee,eps)
IF(numpe==1.AND.i==1)THEN
WRITE(11,'(A,I5)') "Point ",i
; WRITE(11,'(6E12.4)') sigma
END IF
END DO gauss_pts_2
IF (numpe==1) WRITE(11,*)"This analysis took :", elap_time( )-timest(1)
CALL shutdown( )
END PROGRAM p121
Scalar integers:
i
simple counter
element counter
iel
iteration counter
iters
simple counter
j
simple counter
k
iteration ceiling
limit
