FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
solve_LINEQ_direct_serial_lag.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
7 use hecmw_util
11 ! use hecmw_solver_
12
13contains
14
15 subroutine solve_lineq_serial_lag_hecmw_init(hecMAT,fstrMAT,is_sym)
16 implicit none
17 type (hecmwST_matrix) :: hecMAT
18 type (fstrST_matrix_contact_lagrange) :: fstrMAT
19 logical :: is_sym
20
21 call set_pointersandindices_directsolver(hecmat,fstrmat,is_sym)
22
24
25
26 subroutine solve_lineq_serial_lag_hecmw(hecMESH,hecMAT,fstrMAT)
27 implicit none
28 type (hecmwST_local_mesh) :: hecMESH
29 type (hecmwST_matrix) :: hecMAT
30 type (fstrST_matrix_contact_lagrange) :: fstrMAT
31 integer (kind=4) :: ntdf, ilag_sta
32 integer (kind=4) :: numNon0
33 integer (kind=4) :: ierr, nprocs, myrank
34
35 real(kind=8), allocatable :: b(:)
36
37 call hecmw_mat_dump(hecmat, hecmesh)
38
39 call set_values_directsolver(hecmat,fstrmat)
40
41 ntdf = hecmat%NP*hecmat%NDOF + fstrmat%num_lagrange
42 ilag_sta = hecmat%NP*hecmat%NDOF + 1
43 numnon0 = hecmat%NPU*hecmat%NDOF**2+hecmat%NP*hecmat%NDOF*(ntdf+1)/2 &
44 + (fstrmat%numU_lagrange)*hecmat%NDOF+fstrmat%num_lagrange
45
46 allocate(b(size(hecmat%B)))
47 b = hecmat%B
48
49 call hecmw_solve_direct_serial_lag(ntdf, ilag_sta, numnon0, pointers, indices, values, b)
50
51 hecmat%X = b
52
53 deallocate(b)
54
55 call hecmw_mat_dump_solution(hecmat)
56
57 end subroutine solve_lineq_serial_lag_hecmw
58
59
This module provides functions of reconstructing.
subroutine, public hecmw_solve_direct_serial_lag(nrows, ilag_sta, nttbr, pointers, indices, values, b)
I/O and Utility.
Definition: hecmw_util_f.F90:7
This module provides functions to set arrays for direct sparse solver \in the case of using standard ...
subroutine set_values_directsolver(hecmat, fstrmat)
This subroutine sets the array for direct sparse solver that contains \the non-zero items(elements)of...
real(kind=kreal), dimension(:), allocatable values
a
integer(kind=kint), dimension(:), allocatable pointers
ia
subroutine set_pointersandindices_directsolver(hecmat, fstrmat, is_sym)
This subroutine sets index arrays for direct sparse solver from those stored \in the matrix structure...
integer(kind=kint), dimension(:), allocatable indices
ja
subroutine solve_lineq_serial_lag_hecmw(hecmesh, hecmat, fstrmat)
subroutine solve_lineq_serial_lag_hecmw_init(hecmat, fstrmat, is_sym)