FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_EIG_lanczos_util.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
6
7contains
8
10 subroutine lanczos_set_initial_value(hecMESH, hecMAT, fstrEIG, eigvec, p, q, beta)
11 use m_fstr
12 use hecmw_util
13 implicit none
14 type(hecmwst_local_mesh) :: hecMESH
15 type(hecmwst_matrix) :: hecMAT
16 type(fstr_eigen) :: fstrEIG
17 integer(kind=kint) :: N, NP, NDOF, NDOF2, NNDOF, NPNDOF
18 integer(kind=kint) :: i, j
19 real(kind=kreal) :: eigvec(:, :), p(:), beta, chk, sigma
20 real(kind=kreal), allocatable :: temp(:)
21 real(kind=kreal), pointer :: q(:), mass(:), filter(:)
22 logical :: is_free
23
24 n = hecmat%N
25 np = hecmat%NP
26 ndof = hecmesh%n_dof
27 ndof2 = ndof*ndof
28 nndof = n *ndof
29 npndof = np*ndof
30
31 sigma = 0.1d0
32 mass => fstreig%mass
33 filter => fstreig%filter
34
35 allocate(temp(nndof))
36 temp = 0.0d0
37
39 if(fstreig%is_free)then
40 do i = 1, np
41 do j = 1, ndof
42 hecmat%D(ndof2*(i-1) + (ndof+1)*(j-1) + 1) = hecmat%D(ndof2*(i-1) + (ndof+1)*(j-1) + 1) + sigma * mass(ndof*(i-1) + j)
43 enddo
44 enddo
45 endif
46
47 call urand1(nndof, temp, hecmesh%my_rank)
48
49 do i = 1, nndof
50 temp(i) = temp(i) * filter(i)
51 enddo
52
54 do i = 1, nndof
55 eigvec(i,1) = mass(i) * temp(i)
56 enddo
57
58 chk = 0.0d0
59 do i = 1, nndof
60 chk = chk + temp(i) * eigvec(i,1)
61 enddo
62 call hecmw_allreduce_r1(hecmesh, chk, hecmw_sum)
63 beta = dsqrt(chk)
64
65 if(beta == 0.0d0)then
66 call hecmw_finalize()
67 stop "Self-orthogonal"
68 endif
69
70 chk = 1.0d0/beta
71 do i = 1, nndof
72 q(i) = temp(i) * chk
73 enddo
74
75 do i = 1, nndof
76 p(i) = mass(i) * q(i)
77 enddo
78 end subroutine lanczos_set_initial_value
79
81 subroutine evsort(EIG, NEW, NEIG)
82 use hecmw
83 implicit none
84 integer(kind=kint) :: i, j, n, ip, minloc, NEIG, IBAF, NEW(NEIG)
85 real(kind=kreal) :: emin, eig(neig)
86
87 do i = 1, neig
88 new(i) = i
89 enddo
90
91 n = neig-1
92 do i = 1, n
93 minloc = i
94 emin = dabs(eig(new(i)))
95 ip = i+1
96 do j=ip,neig
97 if(dabs(eig(new(j))).LT.emin)then
98 minloc = j
99 emin = dabs(eig(new(j)))
100 endif
101 enddo
102 ibaf = new(i)
103 new(i) = new(minloc)
104 new(minloc) = ibaf
105 enddo
106 end subroutine evsort
107
108 subroutine urand1(N, X, SHIFT)
109 use hecmw
110 implicit none
111 real(kind=kreal) :: x(n), invm
112 integer(kind=kint), parameter :: MM = 1664501
113 integer(kind=kint), parameter :: LAMBDA = 1229
114 integer(kind=kint), parameter :: MU = 351750
115 integer(kind=kint) :: i, N, IR, SHIFT
116
117 ir = 0
118 invm = 1.0d0 / mm
119 do i = 1, shift
120 ir = mod( lambda * ir + mu, mm)
121 enddo
122 do i = shift+1, shift+n
123 ir = mod( lambda * ir + mu, mm)
124 x(i-shift) = invm * ir
125 enddo
126 end subroutine urand1
127
129
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=kint), parameter hecmw_sum
integer(kind=4), parameter kreal
subroutine hecmw_finalize
Definition: hecmw.f90:6
subroutine evsort(eig, new, neig)
Sort eigenvalues.
subroutine urand1(n, x, shift)
subroutine lanczos_set_initial_value(hecmesh, hecmat, fstreig, eigvec, p, q, beta)
Initialize Lanczos iterations.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
Package of data used by Lanczos eigenvalue solver.
Definition: m_fstr.f90:562