28 type (sparse_matrix),
intent(inout) :: spmat
29 type (hecmwst_matrix),
intent(in) :: hecmat
31 type (hecmwst_local_mesh),
intent(in) :: hecmesh
32 integer(kind=kint) :: ndof, ndof2, n_loc, nl, nu, nz, nn,nlag
33 ndof=hecmat%NDOF; ndof2=ndof*ndof
40 n_loc=hecmat%N*ndof+fstrmat%num_lagrange
43 nz=nn*(ndof2+ndof)/2+nu*ndof2 &
44 +fstrmat%numU_lagrange*ndof
49 +(fstrmat%numL_lagrange+fstrmat%numU_lagrange)*ndof
55 if(fstrmat%num_lagrange > 0) &
56 print
'(I3,A,4I10,A,2I10)',
myrank,
' sparse_matrix_init ',hecmat%N,hecmat%NP,n_loc,nz,
' LAG',fstrmat%num_lagrange,spmat%offset
57 nlag = fstrmat%num_lagrange
58 call hecmw_allreduce_i1(hecmesh, nlag, hecmw_sum)
59 if(
myrank == 0) print *,
'total number of contact nodes:',nlag
60 spmat%timelog = hecmat%Iarray(22)
62 call sparse_matrix_para_contact_set_prof(spmat, hecmat, fstrmat)
64 call sparse_matrix_contact_set_prof(spmat, hecmat, fstrmat)
68 subroutine sparse_matrix_contact_set_prof(spMAT, hecMAT, fstrMAT)
69 type(sparse_matrix),
intent(inout) :: spmat
70 type(hecmwst_matrix),
intent(in) :: hecmat
72 integer(kind=kint) :: ndof, ndof2
73 integer(kind=kint) :: m, i, idof, i0, ii, ls, le, l, j, j0, jdof, jdofs
76 ndof=hecmat%NDOF; ndof2=ndof*ndof
81 i0=spmat%offset+ndof*(i-1)
86 ls=hecmat%indexL(i-1)+1
91 j0=spmat%offset+ndof*(j-1)
114 ls=hecmat%indexU(i-1)+1
118 if (j <= hecmat%N)
then
119 j0=spmat%offset+ndof*(j-1)
121 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
133 if (fstrmat%num_lagrange > 0)
then
134 j0=spmat%offset+ndof*hecmat%N
135 ls=fstrmat%indexU_lagrange(i-1)+1
136 le=fstrmat%indexU_lagrange(i)
138 j=fstrmat%itemU_lagrange(l)
147 if (fstrmat%num_lagrange > 0)
then
148 i0=spmat%offset+ndof*hecmat%N
149 do i=1,fstrmat%num_lagrange
153 ls=fstrmat%indexL_lagrange(i-1)+1
154 le=fstrmat%indexL_lagrange(i)
156 j=fstrmat%itemL_lagrange(l)
157 j0=spmat%offset+ndof*(j-1)
173 if (m-1 < spmat%NZ) spmat%NZ=m-1
174 if (m-1 /= spmat%NZ)
then
175 write(*,*)
'ERROR: sparse_matrix_contact_set_prof on rank ',
myrank
176 write(*,*)
'm-1 = ',m-1,
', NZ=',spmat%NZ,
',num_lagrange=',fstrmat%num_lagrange
177 call hecmw_abort(hecmw_comm_get_comm())
179 end subroutine sparse_matrix_contact_set_prof
181 subroutine sparse_matrix_para_contact_set_prof(spMAT, hecMAT, fstrMAT)
182 type(sparse_matrix),
intent(inout) :: spmat
183 type(hecmwst_matrix),
intent(in) :: hecmat
185 integer(kind=kint) :: ndof, ndof2
186 integer(kind=kint) :: m, i, idof, i0, ii, ls, le, l, j, j0, jdof, jdofs
190 write(*,*)
'ERROR: sparse_matrix_para_contact_set_prof on rank ',
myrank
191 write(*,*)
'spMAT%type must be SPARSE_MATRIX_TYPE_COO, only for mumps'
192 call hecmw_abort(hecmw_comm_get_comm())
195 ndof=hecmat%NDOF; ndof2=ndof*ndof
199 if(i <= hecmat%N)
then
201 i0=spmat%offset+ndof*(i-1)
205 ls=hecmat%indexL(i-1)+1
209 j0=spmat%offset+ndof*(j-1)
231 ls=hecmat%indexU(i-1)+1
235 if (j <= hecmat%N)
then
236 j0=spmat%offset+ndof*(j-1)
238 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
251 if (fstrmat%num_lagrange > 0)
then
252 j0=spmat%offset+ndof*hecmat%N
253 ls=fstrmat%indexU_lagrange(i-1)+1
254 le=fstrmat%indexU_lagrange(i)
256 j=fstrmat%itemU_lagrange(l)
266 i0 = spmat%conv_ext(ndof*(i-hecmat%N))-ndof
270 ls=hecmat%indexL(i-1)+1
274 if (j <= hecmat%N)
then
275 j0=spmat%offset+ndof*(j-1)
277 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
301 ls=hecmat%indexU(i-1)+1
305 if (j <= hecmat%N)
then
306 j0=spmat%offset+ndof*(j-1)
308 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
322 if (fstrmat%num_lagrange > 0)
then
323 j0=spmat%offset+ndof*hecmat%N
324 ls=fstrmat%indexU_lagrange(i-1)+1
325 le=fstrmat%indexU_lagrange(i)
330 j=fstrmat%itemU_lagrange(l)
343 if (fstrmat%num_lagrange > 0)
then
344 i0=spmat%offset+ndof*hecmat%N
345 do i=1,fstrmat%num_lagrange
347 ls=fstrmat%indexL_lagrange(i-1)+1
348 le=fstrmat%indexL_lagrange(i)
350 j=fstrmat%itemL_lagrange(l)
351 if (j <= hecmat%N)
then
352 j0=spmat%offset+ndof*(j-1)
354 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
369 if (m-1 /= spmat%NZ)
then
370 write(*,*)
'ERROR: sparse_matrix_para_contact_set_prof on rank ',
myrank
371 write(*,*)
'm-1 = ',m-1,
', NZ=',spmat%NZ,
',num_lagrange=',fstrmat%num_lagrange
372 call hecmw_abort(hecmw_comm_get_comm())
374 end subroutine sparse_matrix_para_contact_set_prof
377 type(sparse_matrix),
intent(inout) :: spmat
378 type(hecmwst_matrix),
intent(in) :: hecmat
380 integer(kind=kint) :: ndof, ndof2
381 integer(kind=kint) :: m, i, idof, i0, ii, ls, le, l, j, j0, jdof, jdofs
382 integer(kind=kint) :: offset_l, offset_d, offset_u
383 ndof=hecmat%NDOF; ndof2=ndof*ndof
388 i0=spmat%offset+ndof*(i-1)
391 stop
"ERROR: sparse_matrix_contact_set_a"
394 ls=hecmat%indexL(i-1)+1
399 j0=spmat%offset+ndof*(j-1)
403 offset_l=ndof2*(l-1)+ndof*(idof-1)
406 if( spmat%IRN(m)/=ii ) stop
"ERROR: sparse_matrix_contact_set_a"
408 if (spmat%JCN(m)/=j0+jdof) stop
"ERROR: sparse_matrix_contact_set_a"
409 spmat%A(m)=hecmat%AL(offset_l+jdof)
415 offset_d=ndof2*(i-1)+ndof*(idof-1)
419 if( spmat%IRN(m)/=ii ) stop
"ERROR: sparse_matrix_contact_set_a"
421 if (spmat%JCN(m)/=i0+jdof) stop
"ERROR: sparse_matrix_contact_set_a"
422 spmat%A(m)=hecmat%D(offset_d+jdof)
426 ls=hecmat%indexU(i-1)+1
430 if (j <= hecmat%N)
then
431 j0=spmat%offset+ndof*(j-1)
433 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
436 offset_u=ndof2*(l-1)+ndof*(idof-1)
439 if( spmat%IRN(m)/=ii ) stop
"ERROR: sparse_matrix_contact_set_a"
441 if (spmat%JCN(m)/=j0+jdof) stop
"ERROR: sparse_matrix_contact_set_a"
442 spmat%A(m)=hecmat%AU(offset_u+jdof)
447 if (fstrmat%num_lagrange > 0)
then
448 j0=spmat%offset+ndof*hecmat%N
449 ls=fstrmat%indexU_lagrange(i-1)+1
450 le=fstrmat%indexU_lagrange(i)
453 if( spmat%IRN(m)/=ii ) stop
"ERROR: sparse_matrix_contact_set_a"
455 if (spmat%JCN(m)/=j0+fstrmat%itemU_lagrange(l)) &
456 stop
"ERROR: sparse_matrix_contact_set_a"
457 spmat%A(m)=fstrmat%AU_lagrange((l-1)*ndof+idof)
464 if (fstrmat%num_lagrange > 0)
then
465 i0=spmat%offset+ndof*hecmat%N
466 do i=1,fstrmat%num_lagrange
470 ls=fstrmat%indexL_lagrange(i-1)+1
471 le=fstrmat%indexL_lagrange(i)
473 j=fstrmat%itemL_lagrange(l)
474 j0=spmat%offset+ndof*(j-1)
477 if( spmat%IRN(m)/=ii ) stop
"ERROR: sparse_matrix_contact_set_a"
479 if (spmat%JCN(m)/=j0+jdof) &
480 stop
"ERROR: sparse_matrix_contact_set_a"
481 spmat%A(m)=fstrmat%AL_lagrange((l-1)*ndof+jdof)
495 stop
"ERROR: sparse_matrix_contact_set_a"
496 if (m-1 /= spmat%NZ) stop
"ERROR: sparse_matrix_contact_set_a"
500 type(sparse_matrix),
intent(inout) :: spmat
501 type(hecmwst_matrix),
intent(in) :: hecmat
503 type(hecmwst_matrix),
intent(in) :: conmat
504 integer(kind=kint) :: ndof, ndof2
505 integer(kind=kint) :: m, i, idof, i0, ii, ls, le, l, j, j0, jdof, jdofs
506 integer(kind=kint) :: offset_l, offset_d, offset_u
509 write(*,*)
'ERROR: sparse_matrix_para_contact_set_vals on rank ',
myrank
510 write(*,*)
'spMAT%type must be SPARSE_MATRIX_TYPE_COO, only for mumps'
511 call hecmw_abort(hecmw_comm_get_comm())
514 ndof=hecmat%NDOF; ndof2=ndof*ndof
518 if(i <= hecmat%N)
then
519 i0=spmat%offset+ndof*(i-1)
524 ls=hecmat%indexL(i-1)+1
528 if (j > hecmat%N) stop
'j > hecMAT%N'
529 j0=spmat%offset+ndof*(j-1)
530 offset_l=ndof2*(l-1)+ndof*(idof-1)
533 stop
"ERROR: sparse_matrix_contact_set_a"
534 if (spmat%JCN(m)/=j0+jdof) &
535 stop
"ERROR: sparse_matrix_contact_set_a"
536 spmat%A(m)=hecmat%AL(offset_l+jdof) + conmat%AL(offset_l+jdof)
543 offset_d=ndof2*(i-1)+ndof*(idof-1)
547 stop
"ERROR: sparse_matrix_contact_set_a"
548 if (spmat%JCN(m)/=i0+jdof) &
549 stop
"ERROR: sparse_matrix_contact_set_a"
551 spmat%A(m)=hecmat%D(offset_d+jdof) + conmat%D(offset_d+jdof)
556 ls=hecmat%indexU(i-1)+1
560 if (j <= hecmat%N)
then
561 j0=spmat%offset+ndof*(j-1)
563 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
566 offset_u=ndof2*(l-1)+ndof*(idof-1)
569 stop
"ERROR: sparse_matrix_contact_set_a"
570 if (spmat%JCN(m)/=j0+jdof) &
571 stop
"ERROR: sparse_matrix_contact_set_a"
572 spmat%A(m)=hecmat%AU(offset_u+jdof) + conmat%AU(offset_u+jdof)
578 if (fstrmat%num_lagrange > 0)
then
579 j0=spmat%offset+ndof*hecmat%N
580 ls=fstrmat%indexU_lagrange(i-1)+1
581 le=fstrmat%indexU_lagrange(i)
584 stop
"ERROR: sparse_matrix_contact_set_a"
585 if (spmat%JCN(m)/=j0+fstrmat%itemU_lagrange(l)) &
586 stop
"ERROR: sparse_matrix_contact_set_a"
587 spmat%A(m)=fstrmat%AU_lagrange((l-1)*ndof+idof)
595 i0 = spmat%conv_ext(ndof*(i-hecmat%N))-ndof
599 ls=hecmat%indexL(i-1)+1
603 if (j <= hecmat%N)
then
604 j0=spmat%offset+ndof*(j-1)
606 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
609 offset_l=ndof2*(l-1)+ndof*(idof-1)
614 stop
"ERROR: sparse_matrix_contact_set_a"
615 if (spmat%JCN(m)/=j0+jdof) &
616 stop
"ERROR: sparse_matrix_contact_set_a"
617 spmat%A(m) = conmat%AL(offset_l+jdof)
624 offset_d=ndof2*(i-1)+ndof*(idof-1)
628 stop
"ERROR: sparse_matrix_contact_set_a"
629 if (spmat%JCN(m)/=i0+jdof) &
630 stop
"ERROR: sparse_matrix_contact_set_a"
633 spmat%A(m) = conmat%D(offset_d+jdof)
639 ls=hecmat%indexU(i-1)+1
643 if (j <= hecmat%N)
then
644 j0=spmat%offset+ndof*(j-1)
646 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
649 offset_u=ndof2*(l-1)+ndof*(idof-1)
652 stop
"ERROR: sparse_matrix_contact_set_a"
653 if (spmat%JCN(m)/=j0+jdof) stop
"ERROR: sparse_matrix_contact_set_a"
654 spmat%A(m) = conmat%AU(offset_u+jdof)
661 if (fstrmat%num_lagrange > 0)
then
662 j0=spmat%offset+ndof*hecmat%N
663 ls=fstrmat%indexU_lagrange(i-1)+1
664 le=fstrmat%indexU_lagrange(i)
670 stop
"ERROR: sparse_matrix_contact_set_a"
671 if (spmat%JCN(m)/=j0+fstrmat%itemU_lagrange(l)) &
672 stop
"ERROR: sparse_matrix_contact_set_a"
673 spmat%A(m)=fstrmat%AU_lagrange((l-1)*ndof+idof)
684 if (fstrmat%num_lagrange > 0)
then
685 i0=spmat%offset+ndof*hecmat%N
686 do i=1,fstrmat%num_lagrange
688 ls=fstrmat%indexL_lagrange(i-1)+1
689 le=fstrmat%indexL_lagrange(i)
691 j=fstrmat%itemL_lagrange(l)
692 if (j <= hecmat%N)
then
693 j0=spmat%offset+ndof*(j-1)
695 j0=spmat%conv_ext(ndof*(j-hecmat%N))-ndof
700 stop
"ERROR: sparse_matrix_contact_set_a"
701 if (spmat%JCN(m)/=j0+jdof) &
702 stop
"ERROR: sparse_matrix_contact_set_a"
703 spmat%A(m)=fstrmat%AL_lagrange((l-1)*ndof+jdof)
712 if(spmat%IRN(ii+1-spmat%offset)/=m) &
713 stop
"ERROR: sparse_matrix_contact_set_a"
715 if (m-1 /= spmat%NZ)
then
716 write(*,*)
'ERROR: sparse_matrix_para_contact_set_vals on rank ',
myrank
717 write(*,*)
'm-1 = ',m-1,
', NZ=',spmat%NZ,
',num_lagrange=',fstrmat%num_lagrange
718 call hecmw_abort(hecmw_comm_get_comm())
724 type (sparse_matrix),
intent(inout) :: spmat
725 type (hecmwst_matrix),
intent(in) :: hecmat
727 integer :: nndof,npndof,ierr,i
728 allocate(spmat%rhs(spmat%N_loc), stat=ierr)
730 write(*,*)
" Allocation error, spMAT%rhs"
731 call hecmw_abort(hecmw_comm_get_comm())
733 nndof=hecmat%N*hecmat%ndof
734 npndof=hecmat%NP*hecmat%ndof
736 spmat%rhs(i) = hecmat%b(i)
739 if (fstrmat%num_lagrange > 0)
then
740 do i=1,fstrmat%num_lagrange
741 spmat%rhs(nndof+i) = hecmat%b(npndof+i)
748 type (sparse_matrix),
intent(inout) :: spmat
749 type (hecmwst_matrix),
intent(in) :: hecmat
751 type (hecmwst_matrix),
intent(in) :: conmat
752 integer :: nndof,npndof,ierr,ndof,i,i0,j
753 real(kind=kreal),
allocatable :: rhs_con(:), rhs_con_sum(:)
756 allocate(rhs_con(spmat%N),stat=ierr)
759 do i= hecmat%N+1,hecmat%NP
760 i0=spmat%conv_ext(ndof*(i-hecmat%N))-ndof
761 if((i0 < 0.or.i0 > spmat%N))
then
763 if(conmat%b((i-1)*ndof+j) /= 0.0d0)
then
764 print *,
myrank,
'i0',i,spmat%N,
'conMAT%b',conmat%b((i-1)*ndof+j)
769 if(i0 > spmat%N - ndof)
then
770 print *,
myrank,
'ext out',hecmat%N,hecmat%NP,i,i0
773 if(conmat%b((i-1)*ndof+j) /= 0.0d0)
then
774 rhs_con(i0+j) = conmat%b((i-1)*ndof+j)
779 allocate(rhs_con_sum(spmat%N))
780 call hecmw_allreduce_dp(rhs_con, rhs_con_sum, spmat%N, hecmw_sum, hecmw_comm_get_comm())
783 allocate(spmat%rhs(spmat%N_loc), stat=ierr)
785 write(*,*)
" Allocation error, spMAT%rhs"
786 call hecmw_abort(hecmw_comm_get_comm())
788 nndof=hecmat%N*hecmat%ndof
789 npndof=hecmat%NP*hecmat%ndof
791 spmat%rhs(i) = rhs_con_sum(spmat%offset+i) + hecmat%b(i) + conmat%b(i)
793 deallocate(rhs_con_sum)
794 if (fstrmat%num_lagrange > 0)
then
795 do i=1,fstrmat%num_lagrange
796 spmat%rhs(nndof+i) = conmat%b(npndof+i)
803 type (sparse_matrix),
intent(inout) :: spmat
804 type (hecmwst_matrix),
intent(inout) :: hecmat
806 integer :: nndof,npndof,i
807 nndof=hecmat%N*hecmat%ndof
808 npndof=hecmat%NP*hecmat%ndof
810 hecmat%x(i) = spmat%rhs(i)
813 if (fstrmat%num_lagrange > 0)
then
814 do i=1,fstrmat%num_lagrange
815 hecmat%x(npndof+i) = spmat%rhs(nndof+i)
818 deallocate(spmat%rhs)
This module defined coomon data and basic structures for analysis.
integer(kind=kint) myrank
PARALLEL EXECUTION.
logical paracontactflag
PARALLEL CONTACT FLAG.
This module provides conversion routines between HEC data structure and DOF based sparse matrix struc...
subroutine, public sparse_matrix_hec_set_conv_ext(spMAT, hecMESH, ndof)
This module provides DOF based sparse matrix data structure (CSR and COO)
subroutine, public sparse_matrix_init(spMAT, N_loc, NZ)
logical function, public sparse_matrix_is_sym(spMAT)
integer(kind=kint), parameter, public sparse_matrix_type_coo
integer(kind=kint), parameter, public sparse_matrix_type_csr