25 integer(kind=kint) :: N, NP, NDOF, NNDOF, NPNDOF
26 integer(kind=kint) :: iter, maxiter, nget, ierr
27 integer(kind=kint) :: i, j, k, in, jn, kn, ik, it
28 integer(kind=kint) :: ig, ig0, is0, ie0, its0, ite0
29 real(kind=
kreal) :: t1, t2, tolerance
30 real(kind=
kreal) :: alpha, beta, beta0
31 real(kind=
kreal),
allocatable :: s(:), t(:), p(:)
32 logical :: is_converge
40 allocate(fstreig%filter(npndof))
41 fstreig%filter = 1.0d0
44 do ig0 = 1, fstrsolid%BOUNDARY_ngrp_tot
45 ig = fstrsolid%BOUNDARY_ngrp_ID(ig0)
46 is0 = hecmesh%node_group%grp_index(ig-1) + 1
47 ie0 = hecmesh%node_group%grp_index(ig )
48 it = fstrsolid%BOUNDARY_ngrp_type(ig0)
49 its0 = (it - mod(it,10))/10
53 in = hecmesh%node_group%grp_item(ik)
54 if(ndof < ite0) ite0 = ndof
57 fstreig%filter((in-1)*ndof+i) = 0.0d0
62 call hecmw_allreduce_i1(hecmesh, jn,
hecmw_sum)
64 fstreig%is_free = .true.
66 write(*,*)
'** free modal analysis: shift factor = 0.1'
70 call hecmw_update_m_r(hecmesh, fstreig%filter, np, ndof)
74 if(fstreig%filter(i) == 1.0d0) in = in + 1
76 call hecmw_allreduce_i1(hecmesh, in,
hecmw_sum)
78 fstreig%maxiter = fstreig%maxiter + 1
79 if(in < fstreig%maxiter)
then
81 write(
imsg,*)
'** changed maxiter to system matrix size.'
86 if(in < fstreig%nget)
then
90 maxiter = fstreig%maxiter
92 allocate(q(0:maxiter))
93 allocate(q(0)%q(npndof))
94 allocate(q(1)%q(npndof))
95 allocate(fstreig%eigval(maxiter))
96 allocate(fstreig%eigvec(npndof, maxiter))
97 allocate(tri%alpha(maxiter))
98 allocate(tri%beta(maxiter))
103 fstreig%eigval = 0.0d0
104 fstreig%eigvec = 0.0d0
116 hecmat%Iarray(98) = 1
117 hecmat%Iarray(97) = 1
121 write(
imsg,*)
' ***** STAGE Begin Lanczos loop **'
124 do iter = 1, maxiter-1
130 call solve_lineq(hecmesh, hecmat)
132 allocate(q(iter+1)%q(npndof))
135 t(i) = hecmat%X(i) * fstreig%filter(i)
141 t(i) = t(i) - tri%beta(iter) * q(iter-1)%q(i)
146 alpha = alpha + p(i) * t(i)
148 call hecmw_allreduce_r1(hecmesh, alpha,
hecmw_sum)
149 tri%alpha(iter) = alpha
153 t(i) = t(i) - tri%alpha(iter) * q(iter)%q(i)
160 s(i) = fstreig%mass(i) * t(i)
166 t1 = t1 + q(j)%q(i) * s(i)
168 call hecmw_allreduce_r1(hecmesh, t1,
hecmw_sum)
170 t(i) = t(i) - t1 * q(j)%q(i)
176 s(i) = fstreig%mass(i) * t(i)
181 beta = beta + s(i) * t(i)
183 call hecmw_allreduce_r1(hecmesh, beta,
hecmw_sum)
184 tri%beta(iter+1) = dsqrt(beta)
188 beta = 1.0d0/tri%beta(iter+1)
191 q(iter+1)%q(i) = t(i) * beta
195 if(iter == 1) beta0 = tri%beta(iter+1)
197 call tridiag(hecmesh, hecmat, fstreig, q, tri, iter, is_converge)
203 if(
associated(q(i)%q))
deallocate(q(i)%q)
205 deallocate(tri%alpha)
215 write(
imsg,*)
' * STAGE Output and postprocessing **'
216 write(
idbg,
'(a,f10.2)')
'Lanczos loop (sec) :', t2 - t1
integer(kind=kint), parameter hecmw_sum
integer(kind=4), parameter kreal
real(kind=kreal) function hecmw_wtime()
This modules just summarizes all modules used in eigen analysis.
subroutine lanczos_set_initial_value(hecMESH, hecMAT, fstrEIG, eigvec, p, q, beta)
Initialize Lanczos iterations.
Lanczos iteration calculation.
subroutine fstr_solve_lanczos(hecMESH, hecMAT, fstrSOLID, fstrEIG)
SOLVE EIGENVALUE PROBLEM.
This module provides a subroutine to find the eigenvalues and eigenvectors of a symmetric tridiagonal...
subroutine tridiag(hecMESH, hecMAT, fstrEIG, Q, Tri, iter, is_converge)
This module defined coomon data and basic structures for analysis.
integer(kind=kint) myrank
PARALLEL EXECUTION.
integer(kind=kint), parameter imsg
integer(kind=kint), parameter idbg
Package of data used by Lanczos eigenvalue solver.