FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
fstr_solve_eigen.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 !-------------------------------------------------------------------------------
7 contains
8 
10  subroutine fstr_solve_eigen( hecMESH, hecMAT, fstrEIG, fstrSOLID, &
11  & fstrRESULT, fstrPARAM, fstrMAT)
12  use hecmw_util
13  use m_fstr
15  use m_fstr_addbc
19  use m_static_lib
22 
23  implicit none
24 
25  type(hecmwst_local_mesh) :: hecMESH
26  type(hecmwst_matrix) :: hecMAT
27  type(fstr_solid) :: fstrSOLID
28  type(hecmwst_result_data) :: fstrRESULT
29  type(fstr_param) :: fstrPARAM
30  type(fstr_eigen) :: fstrEIG
31  type(fstrst_matrix_contact_lagrange) :: fstrMAT
32 
33  type(hecmwst_local_mesh), pointer :: hecMESHmpc
34  type(hecmwst_matrix), pointer :: hecMATmpc
35  real(kind=kreal) :: t1, t2
36 
37  t1 = hecmw_wtime()
38 
39  call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
40 
41  fstrsolid%dunode = 0.0d0
42  call fstr_stiffmatrix(hecmesh, hecmat, fstrsolid, 0.0d0, 0.0d0)
43 
44  call hecmw_mpc_mat_ass(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
45  call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
46  call fstr_addbc(1, hecmesh, hecmatmpc, fstrsolid, fstrparam, fstrmat, 2)
47 
48  call setmass(fstrsolid, hecmesh, hecmat, fstreig)
49  call hecmw_mpc_trans_mass(hecmesh, hecmat, fstreig%mass)
50 
51  call fstr_solve_lanczos(hecmeshmpc, hecmatmpc, fstrsolid, fstreig)
52  call hecmw_mpc_tback_eigvec(hecmesh, hecmat, fstreig%iter, fstreig%eigvec)
53 
54  call fstr_eigen_output(hecmesh, hecmat, fstreig)
55 
56  call fstr_eigen_make_result(hecmesh, hecmat, fstreig, fstrresult)
57 
58  call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
59 
60  t2 = hecmw_wtime()
61 
62  if(myrank == 0)then
63  write(imsg,'("### FSTR_SOLVE_EIGEN FINISHED!")')
64  write(*,'("### FSTR_SOLVE_EIGEN FINISHED!")')
65  endif
66 
67  end subroutine fstr_solve_eigen
68 end module m_fstr_solve_eigen
This module provides functions of reconstructing.
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
real(kind=kreal) function hecmw_wtime()
This module provides a function to deal with prescribed displacement.
Definition: fstr_AddBC.f90:7
subroutine fstr_addbc(cstep, hecMESH, hecMAT, fstrSOLID, fstrPARAM, fstrMAT, iter, conMAT)
Add Essential Boundary Conditions.
Definition: fstr_AddBC.f90:14
Lanczos iteration calculation.
subroutine fstr_solve_lanczos(hecMESH, hecMAT, fstrSOLID, fstrEIG)
SOLVE EIGENVALUE PROBLEM.
subroutine fstr_eigen_make_result(hecMESH, hecMAT, fstrEIG, fstrRESULT)
subroutine fstr_eigen_output(hecMESH, hecMAT, fstrEIG)
Set up lumped mass matrix.
subroutine setmass(fstrSOLID, hecMESH, hecMAT, fstrEIG)
This module provides a function to control eigen analysis.
subroutine fstr_solve_eigen(hecMESH, hecMAT, fstrEIG, fstrSOLID, fstrRESULT, fstrPARAM, fstrMAT)
solve eigenvalue probrem
This module provides function to calcualte tangent stiffness matrix.
subroutine, public fstr_stiffmatrix(hecMESH, hecMAT, fstrSOLID, time, tincr)
接線剛性マトリックスを作成するサブルーチン
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint) myrank
PARALLEL EXECUTION.
Definition: m_fstr.f90:80
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
HECMW to FSTR Mesh Data Converter. Convering Conectivity of Element Type 232, 342 and 352.
This modules just summarizes all modules used in static analysis.
Definition: static_LIB.f90:6
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...
Package of data used by Lanczos eigenvalue solver.
Definition: m_fstr.f90:562
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138