!======================================================================! ! ! ! Software Name : FrontISTR Ver. 3.5 ! ! ! ! Module Name : Static Analysis ! ! ! ! Written by K. Sato (Advancesoft) ! ! X. Yuan (Advancesoft) ! ! ! ! Contact address : IIS,The University of Tokyo, CISS ! ! ! ! "Structural Analysis for Large Scale Assembly" ! ! ! !======================================================================! !C !C*** !> CONSTRUCT the GLOBAL STIFF MATRIX !C*** !C module m_static_mat_ass_main implicit none contains subroutine FSTR_MAT_ASS_MAIN (hecMESH, hecMAT, fstrSOLID) use m_fstr type (hecmwST_matrix) :: hecMAT type (hecmwST_local_mesh) :: hecMESH type (fstr_solid) :: fstrSOLID !** Local variables real(kind=kreal) :: xx(20), yy(20), zz(20), stiffness(20*6, 20*6) integer(kind=kint) :: nodLOCAL(20) integer(kind=kint) :: ndof, itype, iS, iE, ic_type, nn, icel, iiS, j integer(kind=kint) :: flag_shell, flag_solid, mixflag ! Gaku Hashimoto (The University of Tokyo, 2014/03/18) < integer(kind=kint) :: isect, cdsys_ID real(kind=kreal) :: coords(3,3) ! > Gaku Hashimoto (The University of Tokyo, 2014/03/18) !C !C +-------+ !C | INIT. | !C +-------+ !C==== call hecmw_mat_clear(hecMAT) hecMAT%X = 0.d0 ! call HECMW_setloglv(HECMW_LOG_DEBUG) !C !C +-------------------------------+ !C | ELEMENT-by-ELEMENT ASSEMBLING | !C | according to ELEMENT TYPE | !C +-------------------------------+ ndof = hecMAT%NDOF do itype= 1, hecMESH%n_elem_type iS= hecMESH%elem_type_index(itype-1) + 1 iE= hecMESH%elem_type_index(itype ) ic_type= hecMESH%elem_type_item(itype) !C** Ignore link elements if (hecmw_is_etype_link(ic_type)) cycle !C** Set number of nodes nn = hecmw_get_max_node(ic_type) !C element loop do icel= iS, iE !C** node ID iiS= hecMESH%elem_node_index(icel-1) do j=1,nn nodLOCAL(j)= hecMESH%elem_node_item (iiS+j) !C** nodal coordinate xx(j)=hecMESH%node(3*nodLOCAL(j)-2) yy(j)=hecMESH%node(3*nodLOCAL(j)-1) zz(j)=hecMESH%node(3*nodLOCAL(j)) enddo !C** Create local stiffness ! Gaku Hashimoto (The University of Tokyo, 2014/05/17) < isect = hecMESH%section_ID(icel) cdsys_ID = hecMESH%section%sect_orien_ID(isect) if( cdsys_ID > 0 ) call get_coordsys( cdsys_ID, hecMESH, fstrSOLID, coords ) ! call fstr_local_stf_create(hecMESH, ndof, ic_type, icel, xx, yy, zz, fstrSOLID%elements(icel)%gausses, & ! fstrSOLID%elements(icel)%iset, stiffness) call fstr_local_stf_create(hecMESH, ndof, ic_type, icel, xx, yy, zz, fstrSOLID%elements(icel)%gausses, & fstrSOLID%elements(icel)%iset, stiffness, cdsys_ID, coords) ! > Gaku Hashimoto (The University of Tokyo, 2014/05/17) !== CONSTRUCT the GLOBAL MATRIX STARTED call hecmw_mat_ass_elem(hecMAT, nn, nodLOCAL, stiffness) enddo enddo !* for EQUATION ! call hecmw_mat_ass_equation ( hecMESH, hecMAT ) end subroutine FSTR_MAT_ASS_MAIN !> Calculate stiff matrix of current element subroutine FSTR_LOCAL_STF_CREATE(hecMESH, ndof, ic_type, icel, xx, yy, zz, gausses, iset, stiffness, cdsys_ID, coords) use m_fstr use m_static_lib use mMechGauss type (hecmwST_local_mesh) :: hecMESH integer(kind=kint) :: ndof, ic_type, icel, iset real(kind=kreal) :: xx(:), yy(:), zz(:), stiffness(:, :) type( tGaussStatus ), intent(in) :: gausses(:) real(kind=kreal) :: ee, pp, thick,ecoord(3,20), coords(3,3) type( tMaterial ), pointer :: material !** Local variables real(kind=kreal) :: local_stf(1830) integer(kind=kint) :: nn, isect, ihead, mixflag ! Gaku Hashimoto (The University of Tokyo, 2014/03/18) < integer(kind=kint) :: cdsys_ID ! > Gaku Hashimoto (The University of Tokyo, 2014/03/18) nn = hecmw_get_max_node(ic_type) ecoord(1,1:nn) = xx(1:nn) ecoord(2,1:nn) = yy(1:nn) ecoord(3,1:nn) = zz(1:nn) material => gausses(1)%pMaterial ee = material%variables(M_YOUNGS) pp = material%variables(M_POISSON) if ( ic_type==241 .or. ic_type==242 .or. & ic_type==231 .or. ic_type==232 .or. ic_type==2322 ) then thick =1.d0 call STF_C2( ic_type,nn,ecoord(1:2,1:nn),gausses(:),thick,stiffness(1:nn*ndof,1:nn*ndof),iset) else if ( ic_type==301 ) then isect= hecMESH%section_ID(icel) ihead = hecMESH%section%sect_R_index(isect-1) thick = hecMESH%section%sect_R_item(ihead+1) call STF_C1( ic_type,nn,ecoord(:,1:nn),thick,gausses(:),stiffness(1:nn*ndof,1:nn*ndof) ) elseif ( ic_type == 781 ) then !for shell-solid mixed analysis isect= hecMESH%section_ID(icel) ihead = hecMESH%section%sect_R_index(isect-1) thick = hecMESH%section%sect_R_item(ihead+1) nn = 4; ndof = 6; ic_type = 741; mixflag = 1 call STF_Shell_MITC(ic_type, nn, ndof, ecoord(1:3, 1:4), gausses(:), stiffness(1:nn*ndof, 1:nn*ndof), thick, mixflag) ic_type = 781; nn = 8; ndof = 3 return elseif ( ic_type == 761 ) then !for shell-solid mixed analysis isect= hecMESH%section_ID(icel) ihead = hecMESH%section%sect_R_index(isect-1) thick = hecMESH%section%sect_R_item(ihead+1) nn = 3; ndof = 6; ic_type = 731; mixflag = 2 call STF_Shell_MITC(ic_type, nn, ndof, ecoord(1:3, 1:3), gausses(:), stiffness(1:nn*ndof, 1:nn*ndof), thick, mixflag) ic_type = 761; nn = 6; ndof = 3 return else if (ic_type==361) then ! Gaku Hashimoto (The University of Tokyo, 2014/03/18) < ! call STF_C3D8IC( ic_type,nn,ecoord(:,1:nn),gausses(:),stiffness(1:nn*ndof,1:nn*ndof)) call STF_C3D8IC( ic_type,nn,ecoord(:,1:nn),gausses(:),stiffness(1:nn*ndof,1:nn*ndof), cdsys_ID, coords) ! > Gaku Hashimoto (The University of Tokyo, 2014/03/18) else if (ic_type==341 .or. ic_type==351 .or. ic_type==361 .or. & ic_type==342 .or. ic_type==352 .or. ic_type==362 ) then call STF_C3( ic_type,nn,ecoord(:,1:nn),gausses(:),stiffness(1:nn*ndof,1:nn*ndof),1.d0, coords) else if( ( ic_type==741 ) .or. ( ic_type==743 ) .or. ( ic_type==731 ) ) then isect= hecMESH%section_ID(icel) ihead = hecMESH%section%sect_R_index(isect-1) thick = hecMESH%section%sect_R_item(ihead+1) mixflag = 0 call STF_Shell_MITC(ic_type, nn, ndof, ecoord(1:3, 1:nn), gausses(:), stiffness(1:nn*ndof, 1:nn*ndof), thick, mixflag) else if ( ic_type==611) then isect= hecMESH%section_ID(icel) ihead = hecMESH%section%sect_R_index(isect-1) call STF_Beam(ic_type,nn,ecoord,hecMESH%section%sect_R_item(ihead+1:),ee, pp,stiffness(1:nn*ndof,1:nn*ndof)) else if(ic_type== 641) THEN isect= hecMESH%section_ID(icel) ihead = hecMESH%section%sect_R_index(isect-1) call STF_Beam_641( ic_type, nn, ecoord,gausses(:), hecMESH%section%sect_R_item(ihead+1:), stiffness(1:nn*ndof,1:nn*ndof) ) else write(*,*) '###ERROR### : Element type not supported for linear static analysis' write(*,*) ' ic_type = ', ic_type call hecmw_abort(hecmw_comm_get_comm()) endif end subroutine FSTR_LOCAL_STF_CREATE end module m_static_mat_ass_main