FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_solve_heat.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!-------------------------------------------------------------------------------
7contains
8
9 subroutine fstr_solve_heat(hecMESH, hecMAT, fstrRESULT, fstrPARAM, fstrHEAT)
10 use m_fstr
11 use m_heat_init
13 use m_heat_io
14 implicit none
15 integer(kind=kint) :: ISTEP
16 type(hecmwst_local_mesh) :: hecMESH
17 type(hecmwst_matrix) :: hecMAT
18 type(hecmwst_result_data) :: fstrRESULT
19 type(fstr_param) :: fstrPARAM
20 type(fstr_heat) :: fstrHEAT
21
22 integer(kind=kint) :: total_step, restart_step_num
23 real(kind=kreal) :: start_time, end_time
24
25 call heat_init(hecmesh, fstrheat)
26 call heat_init_log(hecmesh)
27
28 if(fstrheat%restart_nout < 0) then
29 call heat_input_restart(fstrheat, hecmesh, restart_step_num, total_step, start_time)
30 end_time = fstrheat%STEP_EETIME(restart_step_num)
31 if( (end_time - start_time) / end_time < 1.d-12 ) then
32 restart_step_num = restart_step_num + 1
33 start_time = 0.0d0
34 endif
35 else
36 restart_step_num = 1
37 total_step = 1
38 start_time = 0.0d0
39 endif
40
41 do istep = restart_step_num, fstrheat%STEPtot
42 fstrheat%is_steady = 0
43 if(fstrheat%STEP_DLTIME(istep) <= 0.0d0) fstrheat%is_steady = 1
44
45 if(hecmesh%my_rank == 0)then
46 write(imsg,"(a,i8,a,i8)")"* Current step / Total step: ", istep, "/", fstrheat%STEPtot
47 write(imsg,"(a,i8)")"* max iteration at each step: ", fstrparam%ITMAX(istep)
48 if(fstrheat%is_steady == 1)then
49 write(imsg,"(a)")"* Steady state analysis"
50 else
51 write(imsg,"(a)")"* Transient anslysis"
52 endif
53 endif
54
55 call heat_solve_tran(hecmesh, hecmat, fstrresult, fstrparam, fstrheat, istep, total_step, start_time)
56 enddo
57
58 call heat_finalize(fstrheat)
59 end subroutine fstr_solve_heat
60end module m_fstr_solve_heat
This module provides a function to control heat analysis.
subroutine fstr_solve_heat(hecmesh, hecmat, fstrresult, fstrparam, fstrheat)
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
This module provides functions to initialize heat analysis.
Definition: heat_init.f90:6
subroutine heat_init_log(hecmesh)
Definition: heat_init.f90:48
subroutine heat_finalize(fstrheat)
Definition: heat_init.f90:63
subroutine heat_init(hecmesh, fstrheat)
Definition: heat_init.f90:10
This module provides a function to control heat analysis.
Definition: heat_io.f90:6
subroutine heat_input_restart(fstrheat, hecmesh, istep, tstep, tt)
Definition: heat_io.f90:11
This module provide a function to control nonsteady heat analysis.
subroutine heat_solve_tran(hecmesh, hecmat, fstrresult, fstrparam, fstrheat, istep, total_step, current_time)
Data for HEAT ANSLYSIS (fstrHEAT)
Definition: m_fstr.f90:394
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138