FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
heat_solve_main.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 heat_solve_main(hecMESH, hecMAT, hecMESHmpc, hecMATmpc, fstrPARAM, fstrHEAT, ISTEP, iterALL, next_time, delta_time)
10 use m_fstr
15 implicit none
16 integer(kind=kint) :: i, iterALL, ISTEP, bup_n_dof
17 real(kind=kreal) :: delta_time, next_time
18 type(hecmwst_local_mesh) :: hecmesh
19 type(hecmwst_matrix) :: hecMAT
20 type(fstr_heat) :: fstrHEAT
21 type(fstr_param) :: fstrPARAM
22 type(hecmwst_local_mesh), pointer :: hecMESHmpc
23 type(hecmwst_matrix), pointer :: hecMATmpc
24 logical :: is_congerged
25
26 iterall = 0
27 do
28 iterall = iterall + 1
29 hecmat%X = 0.0d0
30
31 call heat_mat_ass_conductivity(hecmesh, hecmat, fstrheat, fstrheat%beta)
32 if(fstrheat%is_steady == 0) call heat_mat_ass_capacity(hecmesh, hecmat, fstrheat, delta_time)
33 call heat_mat_ass_boundary(hecmesh, hecmat, hecmeshmpc, hecmatmpc, fstrheat, next_time, delta_time)
34
35 hecmatmpc%Iarray(97) = 1 !Need numerical factorization
36 bup_n_dof = hecmesh%n_dof
37 hecmesh%n_dof = 1
38 call solve_lineq(hecmeshmpc, hecmatmpc)
39 hecmesh%n_dof = bup_n_dof
40 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
41
42 do i = 1, hecmesh%n_node
43 fstrheat%TEMPC(i) = fstrheat%TEMP(i)
44 fstrheat%TEMP (i) = hecmat%X(i)
45 enddo
46
47 call heat_check_convergence(hecmesh, fstrheat, fstrparam, istep, iterall, is_congerged)
48
49 if(is_congerged .or. fstrheat%is_iter_max_limit) exit
50 enddo
51 end subroutine heat_solve_main
52
53 subroutine heat_check_convergence(hecMESH, fstrHEAT, fstrPARAM, ISTEP, iterALL, is_congerged)
54 use m_fstr
55 implicit none
56 integer(kind=kint) :: i, iterALL, iter_max, ISTEP
57 real(kind=kreal) :: val, iter_eps
58 type(hecmwst_local_mesh) :: hecmesh
59 type(fstr_heat) :: fstrHEAT
60 type(fstr_param) :: fstrPARAM
61 logical :: is_congerged
62
63 is_congerged = .false.
64 fstrheat%is_iter_max_limit = .false.
65 iter_max = fstrparam%itmax(istep)
66 iter_eps = fstrparam%eps(istep)
67
68 val = 0.0d0
69 do i = 1, hecmesh%nn_internal
70 val = val + (fstrheat%TEMP(i) - fstrheat%TEMPC(i))**2
71 enddo
72 call hecmw_allreduce_r1(hecmesh, val, hecmw_sum)
73 val = dsqrt(val)
74
75 if(hecmesh%my_rank == 0)then
76 write(*,'(i8,1pe16.6)') iterall, val
77 endif
78
79 if(val < iter_eps)then
80 if(hecmesh%my_rank == 0) write(imsg,*) ' !!! CONVERGENCE ACHIEVED '
81 is_congerged = .true.
82 endif
83
84 if(iter_max <= iterall)then
85 if(hecmesh%my_rank == 0) write(*,*) ' !!! ITERATION COUNT OVER : MAX = ', iter_max
86 fstrheat%is_iter_max_limit = .true.
87 endif
88 end subroutine heat_check_convergence
89end module m_heat_solve_main
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 a subroutine for all boundary conditions needed in heat anaylsis.
subroutine heat_mat_ass_boundary(hecmesh, hecmat, hecmeshmpc, hecmatmpc, fstrheat, next_time, delta_time)
This module provides a subroutine to assemble heat capacity matrix.
subroutine heat_mat_ass_capacity(hecmesh, hecmat, fstrheat, delta_time)
subroutine heat_mat_ass_conductivity(hecmesh, hecmat, fstrheat, beta)
This module provides a function for stationary heat analysis.
subroutine heat_solve_main(hecmesh, hecmat, hecmeshmpc, hecmatmpc, fstrparam, fstrheat, istep, iterall, next_time, delta_time)
subroutine heat_check_convergence(hecmesh, fstrheat, fstrparam, istep, iterall, is_congerged)
subroutine, public solve_lineq(hecmesh, hecmat)
Definition: solve_LINEQ.f90:16
Data for HEAT ANSLYSIS (fstrHEAT)
Definition: m_fstr.f90:394
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138