make_result.f90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
! Copyright (c) 2019 FrontISTR Commons
! This software is released under the MIT License, see LICENSE.txt
!-------------------------------------------------------------------------------
module m_make_result
  private
 
  public:: fstr_write_result
  public:: fstr_write_result_add
  public:: fstr_make_result
  public:: fstr_merge_result
  public:: fstr_reorder_node_shell
  public:: fstr_reorder_rot_shell
  public:: fstr_reorder_node_beam
  public:: setup_contact_output_variables
 
contains
 
  !C***
  !C***
  subroutine fstr_write_result( hecMESH, fstrSOLID, fstrPARAM, istep, time, flag, fstrDYNAMIC)
    use m_fstr
    use hecmw_util
 
    implicit none
    type (hecmwst_local_mesh) :: hecmesh
    type (fstr_solid)         :: fstrsolid
    type (fstr_param       )  :: fstrparam
    integer(kind=kint)        :: istep, flag
    type (fstr_dynamic), intent(in), optional  :: fstrdynamic
    real(kind=kreal)          :: time        
    integer(kind=kint) :: id
    real(kind=kreal)   :: work(1)
    character(len=HECMW_HEADER_LEN) :: header
    character(len=HECMW_MSG_LEN)    :: comment
    character(len=HECMW_NAME_LEN)   :: label, nameid, addfname
 
    ! --- INITIALIZE
    header = '*fstrresult'
    if( present(fstrdynamic) ) then
      comment = 'dynamic_result'
    else
      comment = 'static_result'
    endif
    call hecmw_result_init( hecmesh, istep, header, comment )
 
    ! --- TIME
    id = hecmw_result_dtype_global
    label = 'TOTALTIME'
    work(1) = time
    call hecmw_result_add( id, 1, label, work )
 
    if( present(fstrdynamic) ) then
      call fstr_write_result_add( hecmesh, fstrsolid, istep, "", fstrdynamic )
    else
      call fstr_write_result_add( hecmesh, fstrsolid, istep, "" )
    endif
 
    ! --- WRITE
    nameid = 'fstrRES'
    if( flag==0 ) then
      call hecmw_result_write_by_name( nameid )
    else
      addfname = '_dif'
      call hecmw_result_write_by_addfname( nameid, addfname )
    endif
 
    ! --- FINALIZE
    call hecmw_result_finalize
  end subroutine fstr_write_result
 
  !C***
  !C***
  subroutine fstr_write_result_add( hecMESH, fstrSOLID, istep, label_suffix, fstrDYNAMIC)
    use m_fstr
    use m_out
    use m_static_lib
    use mmaterial
    use hecmw_util
 
    implicit none
    type (hecmwst_local_mesh) :: hecmesh
    type (fstr_solid)         :: fstrsolid
    integer(kind=kint)        :: istep
    character(len=*)          :: label_suffix
    type (fstr_dynamic), intent(in), optional  :: fstrdynamic
    integer(kind=kint) :: n_lyr, ntot_lyr, tmp, is_33shell, is_33beam, cid
    integer(kind=kint) :: i, j, k, ndof, mdof, id, nitem, nn, mm, ngauss, it
    real(kind=kreal), pointer :: tnstrain(:), testrain(:), yield_ratio(:)
    integer(kind=kint) :: idx
    real(kind=kreal), allocatable   :: work(:), unode(:), rnode(:)
    character(len=HECMW_NAME_LEN)   :: s, label, cnum
    character(len=16), allocatable   :: clyr(:)
    logical :: is_dynamic
 
    tnstrain => fstrsolid%TNSTRAIN
    testrain => fstrsolid%TESTRAIN
    yield_ratio => fstrsolid%YIELD_RATIO
 
    is_dynamic = present(fstrdynamic)
 
    if( is_dynamic ) then
      idx = 1
      if( fstrdynamic%idx_eqa==1 .and. istep>0 ) idx = 2
    endif
 
    ndof = hecmesh%n_dof
    mm   = hecmesh%n_node
    if( hecmesh%n_elem > hecmesh%n_node ) mm = hecmesh%n_elem
    if( ndof==2 ) mdof = 3
    if( ndof==3 ) mdof = 6
    if( ndof==4 ) mdof = 6
    if( ndof==6 ) mdof = 6
 
    ntot_lyr   = fstrsolid%max_lyr
    is_33shell = fstrsolid%is_33shell
    is_33beam  = fstrsolid%is_33beam
 
    nn = mm * mdof
    allocate( work(nn) )
 
    ! --- DISPLACEMENT
    if( fstrsolid%output_ctrl(3)%outinfo%on(1) ) then
      if(ndof == 4) then
        id = hecmw_result_dtype_node
        ! for VELOCITY
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(1), 3 )
        allocate( unode(3*hecmesh%n_node) )
        unode = 0.0d0
        do i=1, hecmesh%n_node
          do j = 1, 3
            unode((i-1)*3 + j) = fstrdynamic%DISP((i-1)*4 + j, idx)
          enddo
        enddo
        label = 'VELOCITY'
        call result_add( id, nitem, label, label_suffix, unode )
        deallocate( unode )
        ! for PRESSURE
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(1), 1 )
        allocate( unode(hecmesh%n_node) )
        unode = 0.0d0
        do i=1, hecmesh%n_node
          unode(i) = fstrdynamic%DISP(i*4, idx)
        enddo
        label = 'PRESSURE'
        call result_add( id, nitem, label, label_suffix, unode )
        deallocate( unode )
      else if(ndof == 6) then
        id = hecmw_result_dtype_node
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(1), 3 )
        allocate( unode(hecmesh%n_node*3) )
        unode = 0.0d0
        if( is_dynamic ) then
          do i=1, hecmesh%n_node
            unode((i-1)*3+1:(i-1)*3+3) = fstrdynamic%DISP((i-1)*ndof+1:(i-1)*ndof+3, idx)
          enddo
        else
          do i=1, hecmesh%n_node
            unode((i-1)*3+1:(i-1)*3+3) = fstrsolid%unode((i-1)*ndof+1:(i-1)*ndof+3)
          enddo
        endif
        label = 'DISPLACEMENT'
        call result_add( id, nitem, label, label_suffix, unode )
        deallocate( unode )
      else
        id = hecmw_result_dtype_node
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(1), ndof )
        allocate( unode(hecmesh%n_node*ndof) )
        unode = 0.0d0
        if( is_dynamic ) then
          unode(:) = fstrdynamic%DISP(:,idx)
        else
          unode(:) = fstrsolid%unode
        endif
        label = 'DISPLACEMENT'
        if(is_33beam == 1)then
          call fstr_reorder_node_beam(fstrsolid, hecmesh, unode)
        endif
        if(is_33shell == 1)then
          call fstr_reorder_node_shell(fstrsolid, hecmesh, unode)
        endif
        call result_add( id, nitem, label, label_suffix, unode )
        deallocate( unode )
      endif
    endif
 
    ! --- ROTATION
    if (fstrsolid%output_ctrl(3)%outinfo%on(18)) then
      if(ndof == 6) then
        id = hecmw_result_dtype_node
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(1), 3 )
        label = 'ROTATION'
        allocate( rnode(hecmesh%n_node*3) )
        rnode = 0.0d0
        if( is_dynamic ) then
          do i=1, hecmesh%n_node
            rnode((i-1)*3+1:(i-1)*3+3) = fstrdynamic%DISP((i-1)*ndof+4:(i-1)*ndof+6, idx)
          enddo
        else
          do i=1, hecmesh%n_node
            rnode((i-1)*3+1:(i-1)*3+3) = fstrsolid%unode((i-1)*ndof+4:(i-1)*ndof+6)
          enddo
        endif
        call result_add( id, nitem, label, label_suffix, rnode )
        deallocate( rnode )
      else
        if ( is_33shell == 1) then
          id = hecmw_result_dtype_node
          nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(1), ndof )
          label = 'ROTATION'
          allocate( rnode(hecmesh%n_node*ndof) )
          rnode = 0.0d0
          call fstr_reorder_rot_shell(fstrsolid, hecmesh, rnode)
          call result_add( id, nitem, label, label_suffix, rnode )
          deallocate( rnode )
        end if
      end if
    endif
 
    ! --- VELOCITY
    if( is_dynamic .and. fstrsolid%output_ctrl(3)%outinfo%on(15) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(15), ndof )
      label = 'VELOCITY'
      call result_add( id, nitem, label, label_suffix, fstrdynamic%VEL(:,idx) )
    endif
 
    ! --- ACCELERATION
    if( is_dynamic .and. fstrsolid%output_ctrl(3)%outinfo%on(16) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(16), ndof )
      label = 'ACCELERATION'
      call result_add( id, nitem, label, label_suffix, fstrdynamic%ACC(:,idx) )
    endif
 
    ! --- REACTION FORCE
    if( fstrsolid%output_ctrl(3)%outinfo%on(2) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(2), ndof )
      label = 'REACTION_FORCE'
      call result_add( id, nitem, label, label_suffix, fstrsolid%REACTION )
    endif
 
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 
    if(is_33shell == 1 .or. ndof == 6)then
      call fstr_write_result_main( hecmesh, fstrsolid, fstrsolid%SHELL, "                ", label_suffix )
    else
      call fstr_write_result_main( hecmesh, fstrsolid, fstrsolid%SOLID, "                ", label_suffix )
    endif
 
    !laminated shell
    if( associated(fstrsolid%SHELL) .and. fstrsolid%output_ctrl(3)%outinfo%on(27) ) then
      allocate(clyr(2*ntot_lyr))
      do i=1,ntot_lyr
        write(cnum,"(i0)")i
        clyr(2*i-1)="_L"//trim(cnum)//"+"
        clyr(2*i  )="_L"//trim(cnum)//"-"
      enddo
      do i=1,ntot_lyr
        call fstr_write_result_main( hecmesh, fstrsolid, fstrsolid%SHELL%LAYER(i)%PLUS,  clyr(2*i-1), label_suffix )
        call fstr_write_result_main( hecmesh, fstrsolid, fstrsolid%SHELL%LAYER(i)%MINUS, clyr(2*i  ), label_suffix )
      enddo
      deallocate(clyr)
    endif
 
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    ! --- STRAIN @gauss
    if( fstrsolid%output_ctrl(3)%outinfo%on(9) .and. ndof/=6 ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(9), ndof )
      ngauss = fstrsolid%maxn_gauss
      work(:) = 0.d0
      do k = 1, ngauss
        write(s,*) k
        write(label,'(a,a)') 'GaussSTRAIN',trim(adjustl(s))
        label = adjustl(label)
        do i = 1, hecmesh%n_elem
          if( associated(fstrsolid%elements(i)%gausses) ) then
            if( k <= size(fstrsolid%elements(i)%gausses) ) then
              do j = 1, nitem
                work(nitem*(i-1)+j) = fstrsolid%elements(i)%gausses(k)%strain_out(j)
              enddo
            endif
          end if
        enddo
        call result_add( id, nitem, label, label_suffix, work )
      enddo
    endif
 
    ! --- STRESS @gauss
    if( fstrsolid%output_ctrl(3)%outinfo%on(10) .and. ndof/=6 ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(10), ndof )
      ngauss = fstrsolid%maxn_gauss
      work(:) = 0.d0
      do k = 1, ngauss
        write(s,*) k
        write(label,'(a,a)') 'GaussSTRESS',trim(adjustl(s))
        label = adjustl(label)
        do i = 1, hecmesh%n_elem
          if( associated(fstrsolid%elements(i)%gausses) ) then
            if( k <= size(fstrsolid%elements(i)%gausses) ) then
              do j = 1, nitem
                work(nitem*(i-1)+j) = fstrsolid%elements(i)%gausses(k)%stress_out(j)
              enddo
            endif
          end if
        enddo
        call result_add( id, nitem, label, label_suffix, work )
      enddo
    endif
 
    ! --- PLASTIC STRAIN @gauss
    if( fstrsolid%output_ctrl(3)%outinfo%on(11) .and. fstrsolid%StaticType/=3 ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(11), ndof )
      ngauss = fstrsolid%maxn_gauss
      work(:) = 0.d0
      do k = 1, ngauss
        write(s,*) k
        write(label,'(a,a)') 'PLASTIC_GaussSTRAIN',trim(adjustl(s))
        label = adjustl(label)
        do i = 1, hecmesh%n_elem
          if( associated(fstrsolid%elements(i)%gausses) ) then
            if( k <= size(fstrsolid%elements(i)%gausses) ) then
            work(i) = fstrsolid%elements(i)%gausses(k)%plstrain
          endif
          endif
        enddo
        call result_add( id, nitem, label, label_suffix, work )
      enddo
    endif
 
    ! --- THERMAL STRAIN @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(12) .and. associated(tnstrain) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(12), ndof )
      label = 'THERMAL_NodalSTRAIN'
      call result_add( id, nitem, label, label_suffix, tnstrain )
    endif
 
    ! --- THERMAL STRAIN @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(13) .and. associated(testrain) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(13), ndof )
      label = 'THERMAL_ElementalSTRAIN'
      call result_add( id, nitem, label, label_suffix, testrain )
    endif
 
    ! --- THERMAL STRAIN @gauss
    if( fstrsolid%output_ctrl(3)%outinfo%on(14) .and. associated(testrain) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(14), ndof )
      ngauss = fstrsolid%maxn_gauss
      do k = 1, ngauss
        write(s,*) k
        write(label,'(a,a)') 'THERMAL_GaussSTRAIN',trim(adjustl(s))
        label = adjustl(label)
        do i = 1, hecmesh%n_elem
          if( k > ngauss ) then
            do j = 1, nitem
              work(nitem*(i-1)+j) = 0.d0
            enddo
          else
            do j = 1, nitem
              !                work(nitem*(i-1)+j) = fstrSOLID%elements(i)%gausses(k)%tstrain(j)
            enddo
          end if
        enddo
        call result_add( id, nitem, label, label_suffix, work )
      enddo
    endif
 
    ! --- YIELD RATIO
    if( fstrsolid%output_ctrl(3)%outinfo%on(29) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(29), ndof )
      label = "YIELD_RATIO"
      call result_add( id, nitem, label, label_suffix, yield_ratio )
    endif
 
    ! --- CONTACT NORMAL FORCE @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(30) .and. associated(fstrsolid%CONT_NFORCE) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(30), ndof )
      label = 'CONTACT_NFORCE'
      call result_add( id, nitem, label, label_suffix, fstrsolid%CONT_NFORCE )
    endif
 
    ! --- CONTACT FRICTION FORCE @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(31) .and. associated(fstrsolid%CONT_FRIC) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(31), ndof )
      label = 'CONTACT_FRICTION'
      call result_add( id, nitem, label, label_suffix, fstrsolid%CONT_FRIC )
    endif
 
    ! --- CONTACT RELATIVE VELOCITY @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(32) .and. associated(fstrsolid%CONT_RELVEL) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(32), ndof )
      label = 'CONTACT_RELVEL'
      call result_add( id, nitem, label, label_suffix, fstrsolid%CONT_RELVEL )
    endif
 
    ! --- CONTACT STATE @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(33) .and. associated(fstrsolid%CONT_STATE) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(33), ndof )
      label = 'CONTACT_STATE'
      call result_add( id, nitem, label, label_suffix, fstrsolid%CONT_STATE )
    endif
 
    ! --- CONTACT NORMAL TRACTION @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(36) .and. associated(fstrsolid%CONT_NTRAC) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(36), ndof )
      label = 'CONTACT_NTRACTION'
      call result_add( id, nitem, label, label_suffix, fstrsolid%CONT_NTRAC )
    endif
 
    ! --- CONTACT FRICTION TRACTION @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(37) .and. associated(fstrsolid%CONT_FTRAC) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(37), ndof )
      label = 'CONTACT_FTRACTION'
      call result_add( id, nitem, label, label_suffix, fstrsolid%CONT_FTRAC )
    endif
 
    ! --- TEMPERATURE @node
    !if(present(fstrPARAM))then
    !  if( fstrSOLID%output_ctrl(4)%outinfo%on(38) .and. fstrPARAM%solution_type == kstHEATSTATIC)then
    !    ncomp = ncomp + 1
    !    nn = n_comp_valtype( fstrSOLID%output_ctrl(4)%outinfo%vtype(38), ndof )
    !    fstrRESULT%nn_dof(ncomp) = nn
    !    fstrRESULT%node_label(ncomp) = 'TEMPERATURE'
    !    do i = 1, hecMESH%n_node
    !      do j = 1, nn
    !        fstrRESULT%node_val_item(nitem*(i-1)+j+iitem) = fstrSOLID%temperature(nn*(i-1)+j)
    !      enddo
    !    enddo
    !    iitem = iitem + nn
    !  endif
    !endif
 
    ! --- ELEMACT flag @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(44) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(44), ndof )
      label = 'ELEMACT'
      work(:) = 0.d0
      do i = 1, hecmesh%n_elem
        if( fstrsolid%elements(i)%elemact_flag /= kelact_inactive ) work(i) = 1.d0
      enddo
      call result_add( id, nitem, label, label_suffix, work )
    endif
 
    deallocate( work )
  end subroutine fstr_write_result_add
 
  subroutine result_add( dtype, n_dof, label, label_suffix, data )
    use m_fstr
    use hecmw_util
 
    implicit none
    integer(kind=kint), intent(in) :: dtype
    integer(kind=kint), intent(in) :: n_dof
    character(len=*), intent(in)   :: label
    character(len=*), intent(in)   :: label_suffix
    real(kind=kreal)               :: data(:)
    character(len=HECMW_NAME_LEN)  :: name
 
    name = trim(label)//trim(label_suffix)
    call hecmw_result_add( dtype, n_dof, name, data )
  end subroutine result_add
 
  subroutine fstr_write_result_main( hecMESH, fstrSOLID, RES, clyr, label_suffix )
    use m_fstr
    use m_out
    use m_static_lib
    use mmaterial
    use hecmw_util
 
    implicit none
    type (hecmwST_local_mesh) :: hecMESH
    type (fstr_solid)         :: fstrSOLID
    type (fstr_solid_physic_val) :: RES
    integer(kind=kint)        :: istep, flag
    integer(kind=kint)        :: n_lyr, cid
 
    character(len=HECMW_HEADER_LEN) :: header
    character(len=HECMW_NAME_LEN)   :: s, label, nameID, addfname
    character(len=16)                :: clyr
    character(len=*)                 :: label_suffix
    character(len=12)                :: cnum
    integer(kind=kint) :: i, j, k, ndof, mdof, id, nitem, nn, mm, ngauss, it
    real(kind=kreal), allocatable   :: work(:)
 
    ndof = hecmesh%n_dof
    allocate( work(hecmesh%n_elem) )
 
    ! --- STRAIN @node
    if (fstrsolid%output_ctrl(3)%outinfo%on(3)) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(3), ndof )
      label = 'NodalSTRAIN'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%STRAIN )
    endif
 
    ! --- STRESS @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(4) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(4), ndof )
      label = 'NodalSTRESS'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%STRESS )
    endif
 
    ! --- MISES @node
    if( fstrsolid%output_ctrl(3)%outinfo%on(5) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(5), ndof )
      label = 'NodalMISES'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%MISES )
    endif
 
    ! --- NODAL PRINC STRESS
    if( fstrsolid%output_ctrl(3)%outinfo%on(19) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(19), ndof )
      label = 'NodalPrincipalSTRESS'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%PSTRESS )
    endif
 
    ! --- NODAL PRINC STRAIN
    if( fstrsolid%output_ctrl(3)%outinfo%on(21) ) then
      id = hecmw_result_dtype_node
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(21), ndof )
      label = 'NodalPrincipalSTRAIN'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%PSTRAIN )
    endif
 
    ! --- NODAL PRINC STRESS VECTOR
    if( fstrsolid%output_ctrl(3)%outinfo%on(23) ) then
      id = hecmw_result_dtype_node
      do k=1,3
        write(cnum,'(i0)')k
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(23), ndof )
        label = 'NodalPrincipalSTRESSVector'//trim(cnum)//trim(clyr)
        call result_add( id, nitem, label, label_suffix, res%PSTRESS_VECT(:,k) )
      end do
    endif
 
    ! --- NODAL PRINC STRAIN VECTOR
    if( fstrsolid%output_ctrl(3)%outinfo%on(25) ) then
      id = hecmw_result_dtype_node
      do k=1,3
        write(cnum,'(i0)')k
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(25), ndof )
        label = 'NodalPrincipalSTRAINVector'//trim(cnum)//trim(clyr)
        call result_add( id, nitem, label, label_suffix, res%PSTRAIN_VECT(:,k) )
      end do
    endif
 
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    ! --- STRAIN @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(6) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(6), ndof )
      label = 'ElementalSTRAIN'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%ESTRAIN )
    endif
 
    ! --- STRESS @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(7) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(7), ndof )
      label = 'ElementalSTRESS'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%ESTRESS )
    endif
 
    ! --- NQM @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(35) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(35), ndof )
      label = 'ElementalNQM'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%ENQM )
    endif
 
    ! --- MISES @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(8)) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(8), ndof )
      label = 'ElementalMISES'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%EMISES )
    endif
 
    ! --- Principal_STRESS @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(20) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(20), ndof )
      label = 'ElementalPrincipalSTRESS'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%EPSTRESS )
    endif
 
    ! --- Principal_STRAIN @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(22) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(22), ndof )
      label = 'ElementalPrincipalSTRAIN'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%EPSTRAIN )
    endif
 
    ! --- ELEM PRINC STRESS VECTOR
    if( fstrsolid%output_ctrl(3)%outinfo%on(24) ) then
      id = hecmw_result_dtype_elem
      do k=1,3
        write(cnum,'(i0)')k
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(24), ndof )
        label = 'ElementalPrincipalSTRESSVector'//trim(cnum)//trim(clyr)
        call result_add( id, nitem, label, label_suffix, res%EPSTRESS_VECT(:,k) )
      end do
    endif
 
    !ELEM PRINC STRAIN VECTOR
    if( fstrsolid%output_ctrl(3)%outinfo%on(26) ) then
      id = hecmw_result_dtype_elem
      do k=1,3
        write(cnum,'(i0)')k
        nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(26), ndof )
        label = 'ElementalPrincipalSTRAINVector'//trim(cnum)//trim(clyr)
        call result_add( id, nitem, label, label_suffix, res%EPSTRAIN_VECT(:,k) )
      end do
    endif
 
    ! --- PLASTIC STRAIN @element
    if( fstrsolid%output_ctrl(3)%outinfo%on(43) ) then
      id = hecmw_result_dtype_elem
      nitem = n_comp_valtype( fstrsolid%output_ctrl(3)%outinfo%vtype(43), ndof )
      label = 'ElementalPLSTRAIN'//trim(clyr)
      call result_add( id, nitem, label, label_suffix, res%EPLSTRAIN )
    endif
    deallocate( work )
  
  end subroutine fstr_write_result_main
 
  !C***
  !C***
  subroutine fstr_make_result( hecMESH, fstrSOLID, fstrRESULT, istep, time, fstrDYNAMIC, label_suffix )
    use m_fstr
    use hecmw_util
 
    implicit none
    character(len=*), intent(in), optional  :: label_suffix
    type (hecmwst_local_mesh) :: hecmesh
    type (fstr_solid)         :: fstrsolid
    type(hecmwst_result_data) :: fstrresult
    integer(kind=kint)        :: istep
    real(kind=kreal) :: time, coords(3,3)
    type(fstr_dynamic), intent(in), optional  :: fstrdynamic
    integer(kind=kint) :: n_lyr, ntot_lyr, it, coef33, is_33shell, is_33beam
    integer(kind=kint) :: i, j, k, ndof, mdof, gcomp, gitem, ncomp, nitem, iitem, ecomp, eitem, jitem, nn, mm
    integer(kind=kint) :: idx
    real(kind=kreal), pointer :: tnstrain(:), testrain(:)
    real(kind=kreal), allocatable   ::unode(:)
    character(len=12) :: cnum
    character(len=16), allocatable   :: clyr(:)
    logical :: is_dynamic
 
    is_dynamic = present(fstrdynamic)
 
    tnstrain => fstrsolid%TNSTRAIN
    testrain => fstrsolid%TESTRAIN
 
    ntot_lyr   = fstrsolid%max_lyr
    is_33shell = fstrsolid%is_33shell
    is_33beam  = fstrsolid%is_33beam
 
    mm = hecmesh%n_node
    if( hecmesh%n_elem>hecmesh%n_node ) mm = hecmesh%n_elem
 
    if( is_dynamic ) then
      idx = 1
      if( fstrdynamic%idx_eqa==1 .and. istep>0 ) idx = 2
    endif
 
    ndof = hecmesh%n_dof
    if( ndof==2 ) mdof = 3
    if( ndof==3 ) mdof = 6
    if( ndof==4 ) mdof = 6
    if( ndof==6 ) mdof = 6
 
    if(is_33shell == 1 .and. fstrsolid%output_ctrl(4)%outinfo%on(27) )then
      coef33 = 1 + 2*ntot_lyr
    else
      coef33 = 1
    endif
 
    call hecmw_nullify_result_data( fstrresult )
    gcomp = 0
    gitem = 0
    ncomp = 0
    nitem = 0
    ecomp = 0
    eitem = 0
 
    ! --- COUNT SUM OF ALL NITEM
    ! --- TIME
    gcomp = gcomp + 1
    gitem = gitem + 1
    ! --- DISPLACEMENT
    if( fstrsolid%output_ctrl(4)%outinfo%on(1) ) then
      if(ndof == 4) then
        ncomp = ncomp + 1
        nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), 3 )
        ncomp = ncomp + 1
        nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), 1 )
      else if(ndof == 6) then
        ncomp = ncomp + 1
        nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), 3 )
      else
        ncomp = ncomp + 1
        nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), ndof )
      endif
    endif
    ! --- VELOCITY
    if( is_dynamic .and. fstrsolid%output_ctrl(4)%outinfo%on(15) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(15), ndof )
    endif
    ! --- ACCELERATION
    if( is_dynamic .and. fstrsolid%output_ctrl(4)%outinfo%on(16) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(16), ndof )
    endif
    ! --- TEMPERATURE @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(17) .and. associated(fstrsolid%temperature) ) then
     ncomp = ncomp + 1
     nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(17), ndof )
    endif
    ! --- ROTATION (Only for 781 shell)
    if( fstrsolid%output_ctrl(4)%outinfo%on(18) ) then
      if(ndof == 6) then
        ncomp = ncomp + 1
        nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(18), 3 )
      else
        if( is_33shell == 1 ) then
          ncomp = ncomp + 1
          nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(18), ndof )
        endif
      endif
    endif
    ! --- REACTION FORCE
    if( fstrsolid%output_ctrl(4)%outinfo%on(2) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(2), ndof )
    endif
    ! --- STRAIN @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(3) ) then
      ncomp = ncomp + 1*coef33
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(3), ndof )*coef33
    endif
    ! --- STRESS @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(4) ) then
      ncomp = ncomp + 1*coef33
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(4), ndof )*coef33
    endif
    ! --- MISES @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(5) ) then
      ncomp = ncomp + 1*coef33
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(5), ndof )*coef33
    endif
    ! --- Principal Stress @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(19) ) then
      ncomp = ncomp + 1*coef33
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(19), ndof )*coef33
    endif
    ! --- Principal Strain @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(21) ) then
      ncomp = ncomp + 1*coef33
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(21), ndof )*coef33
    endif
    ! --- Principal Stress Vector @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(23) ) then
      ncomp = ncomp + 3*coef33
      nitem = nitem + 3*n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(23), ndof )*coef33
    endif
    ! --- Principal Strain Vector @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(25) ) then
      ncomp = ncomp + 3*coef33
      nitem = nitem + 3*n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(25), ndof )*coef33
    endif
    ! --- THERMAL STRAIN @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(12) .and. associated(tnstrain) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(12), ndof )
    endif
    ! --- CONTACT NORMAL FORCE @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(30) .and. associated(fstrsolid%CONT_NFORCE) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(30), ndof )
    endif
    ! --- CONTACT FRICTION FORCE @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(31) .and. associated(fstrsolid%CONT_FRIC) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(31), ndof )
    endif
    ! --- CONTACT RELATIVE VELOCITY @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(32) .and. associated(fstrsolid%CONT_RELVEL) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(32), ndof )
    endif
    ! --- CONTACT STATE @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(33) .and. associated(fstrsolid%CONT_STATE) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(33), ndof )
    endif
    ! --- CONTACT NORMAL TRACTION @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(36) .and. associated(fstrsolid%CONT_NTRAC) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(36), ndof )
    endif
    ! --- CONTACT FRICTION TRACTION @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(37) .and. associated(fstrsolid%CONT_FTRAC) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(37), ndof )
    endif
    ! --- NODE ID @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(38) ) then
      ncomp = ncomp + 1
      nitem = nitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(38), ndof )
    endif
    ! --- TEMPERATURE @node
    !if( fstrSOLID%output_ctrl(4)%outinfo%on(41) .and. associated(fstrSOLID%CONT_FTRAC) ) then
    !  ncomp = ncomp + 1
    !  nitem = nitem + n_comp_valtype( fstrSOLID%output_ctrl(4)%outinfo%vtype(41), ndof )
    !endif
 
    ! --- STRAIN @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(6) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(6), ndof )
    endif
    ! --- STRESS @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(7) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(7), ndof )
    endif
    ! --- MISES @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(8) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(8), ndof )
    endif
    ! --- Principal Stress @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(20) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(20), ndof )
    endif
    ! --- Principal Strain @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(22) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(22), ndof )
    endif
    ! --- Principal Stress Vector @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(24) ) then
      ecomp = ecomp + 3
      eitem = eitem + 3*n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(24), ndof )
    endif
    ! --- Principal Strain Vector @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(26) ) then
      ecomp = ecomp + 3
      eitem = eitem + 3*n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(26), ndof )
    endif
    ! --- PLASTIC STRAIN @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(43) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(43), ndof )
    endif
    ! --- MATERIAL @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(34) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(34), ndof )
    endif
    ! --- ORIENTATION @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(42) ) then
      ecomp = ecomp + 3*coef33
      eitem = eitem + 3*n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(42), ndof )*coef33
    endif
    ! --- ELEM ID @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(39) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(39), ndof )
    endif
    ! --- SECTION ID @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(40) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(40), ndof )
    endif
    ! --- ELEMACT flag @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(44) ) then
      ecomp = ecomp + 1
      eitem = eitem + n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(44), ndof )
    endif
 
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    fstrresult%ng_component = gcomp
    fstrresult%nn_component = ncomp
    fstrresult%ne_component = ecomp
    allocate( fstrresult%ng_dof(gcomp) )
    allocate( fstrresult%global_label(gcomp) )
    allocate( fstrresult%global_val_item(gitem) )
    allocate( fstrresult%nn_dof(ncomp) )
    allocate( fstrresult%node_label(ncomp) )
    allocate( fstrresult%node_val_item(nitem*hecmesh%n_node) )
    allocate( fstrresult%ne_dof(ecomp) )
    allocate( fstrresult%elem_label(ecomp) )
    allocate( fstrresult%elem_val_item(eitem*hecmesh%n_elem) )
    ncomp = 0
    iitem = 0
    ecomp = 0
    jitem = 0
 
    ! --- TIME
    fstrresult%ng_dof(1) = 1
    fstrresult%global_label(1) = "TOTALTIME"
    fstrresult%global_val_item(1) = time
 
    ! --- DISPLACEMENT
    if (fstrsolid%output_ctrl(4)%outinfo%on(1) ) then
      if(ndof == 4) then
        ! DIPLACEMENT
        ncomp = ncomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), 3 )
        fstrresult%nn_dof(ncomp) = nn
        fstrresult%node_label(ncomp) = 'VELOCITY'
        do i = 1, hecmesh%n_node
          do j = 1, 3
            fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrdynamic%DISP(4*(i-1)+j,idx)
          enddo
        enddo
        iitem = iitem + nn
        ! PRESSURE
        ncomp = ncomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), 1 )
        fstrresult%nn_dof(ncomp) = nn
        fstrresult%node_label(ncomp) = 'PRESSURE'
        do i = 1, hecmesh%n_node
          fstrresult%node_val_item(nitem*(i-1)+1+iitem) = fstrdynamic%DISP(4*i,idx)
        enddo
        iitem = iitem + nn
      else if(ndof == 6) then
        ncomp = ncomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), 3 )
        fstrresult%nn_dof(ncomp) = nn
        fstrresult%node_label(ncomp) = 'DISPLACEMENT'
        allocate( unode(3*hecmesh%n_node) )
        unode = 0.0d0
        if( is_dynamic ) then
          do i=1, hecmesh%n_node
            unode((i-1)*3+1:(i-1)*3+3) = fstrdynamic%DISP((i-1)*ndof+1:(i-1)*ndof+3, idx)
          enddo
        else
          do i=1, hecmesh%n_node
            unode((i-1)*3+1:(i-1)*3+3) = fstrsolid%unode((i-1)*ndof+1:(i-1)*ndof+3)
          enddo
        endif
        do i = 1, hecmesh%n_node
          do j = 1, nn
            fstrresult%node_val_item(nitem*(i-1)+j+iitem) = unode(nn*(i-1)+j)
          enddo
        enddo
        deallocate( unode )
        iitem = iitem + nn
 
      else
        ncomp = ncomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), ndof )
        fstrresult%nn_dof(ncomp) = nn
        fstrresult%node_label(ncomp) = 'DISPLACEMENT'
        allocate( unode(ndof*hecmesh%n_node) )
        unode = 0.0d0
        if( is_dynamic ) then
          unode(:) = fstrdynamic%DISP(:,idx)
        else
          unode(:) = fstrsolid%unode(:)
        endif
        if(is_33beam == 1)then
          call fstr_reorder_node_beam(fstrsolid, hecmesh, unode)
        endif
        if(is_33shell == 1)then
          call fstr_reorder_node_shell(fstrsolid, hecmesh, unode)
        endif
        do i = 1, hecmesh%n_node
          do j = 1, nn
            fstrresult%node_val_item(nitem*(i-1)+j+iitem) = unode(nn*(i-1)+j)
          enddo
        enddo
        deallocate( unode )
        iitem = iitem + nn
      endif
    endif
 
    ! --- VELOCITY
    if( is_dynamic .and. fstrsolid%output_ctrl(4)%outinfo%on(15) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(15), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'VELOCITY'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrdynamic%VEL(nn*(i-1)+j,idx)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- ACCELERATION
    if( is_dynamic .and. fstrsolid%output_ctrl(4)%outinfo%on(16) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(16), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'ACCELERATION'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrdynamic%ACC(nn*(i-1)+j,idx)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- TEMPERATURE
    if( fstrsolid%output_ctrl(4)%outinfo%on(17) .and. associated(fstrsolid%temperature))then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(17), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'TEMPERATURE'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%temperature(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- ROTATION
    if( fstrsolid%output_ctrl(4)%outinfo%on(18) ) then
 
      if(ndof == 6) then
        ncomp = ncomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), 3 )
        fstrresult%nn_dof(ncomp) = nn
        fstrresult%node_label(ncomp) = 'ROTATION'
        allocate( unode(3*hecmesh%n_node) )
        unode = 0.0d0
        if( is_dynamic ) then
          do i=1, hecmesh%n_node
            unode((i-1)*3+1:(i-1)*3+3) = fstrdynamic%DISP((i-1)*ndof+4:(i-1)*ndof+6, idx)
          enddo
        else
          do i=1, hecmesh%n_node
            unode((i-1)*3+1:(i-1)*3+3) = fstrsolid%unode((i-1)*ndof+4:(i-1)*ndof+6)
          enddo
        endif
        do i = 1, hecmesh%n_node
          do j = 1, nn
            fstrresult%node_val_item(nitem*(i-1)+j+iitem) = unode(nn*(i-1)+j)
          enddo
        enddo
        deallocate( unode )
        iitem = iitem + nn
      else
        if ( is_33shell == 1) then
          ncomp = ncomp + 1
          nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(1), ndof )
          fstrresult%nn_dof(ncomp) = nn
          fstrresult%node_label(ncomp) = 'ROTATION'
          allocate( unode(ndof*hecmesh%n_node) )
          unode = 0.0d0
          call fstr_reorder_rot_shell(fstrsolid, hecmesh, unode)
          do i = 1, hecmesh%n_node
            do j = 1, nn
              fstrresult%node_val_item(nitem*(i-1)+j+iitem) = unode(nn*(i-1)+j)
            enddo
          enddo
          deallocate( unode )
          iitem = iitem + nn
        end if
      end if
    endif
 
    ! --- REACTION FORCE
    if( fstrsolid%output_ctrl(4)%outinfo%on(2) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(2), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'REACTION_FORCE'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%REACTION(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    if(is_33shell == 1 .or. ndof == 6)then
      call fstr_make_result_main( hecmesh, fstrsolid, fstrresult, &
        & fstrsolid%SHELL, nitem, iitem, ncomp, eitem, jitem, ecomp, 1, "                " )
    else
      call fstr_make_result_main( hecmesh, fstrsolid, fstrresult, &
        & fstrsolid%SOLID, nitem, iitem, ncomp, eitem, jitem, ecomp, 1, "                " )
    endif
 
    !laminated shell
    if( associated(fstrsolid%SHELL) .and. fstrsolid%output_ctrl(4)%outinfo%on(27) .and. is_33shell == 1 ) then
      allocate(clyr(2*ntot_lyr))
      do i=1,ntot_lyr
        write(cnum,"(i0)")i
        clyr(2*i-1)="_L"//trim(cnum)//"+"
        clyr(2*i  )="_L"//trim(cnum)//"-"
      enddo
      do i=1,ntot_lyr
        call fstr_make_result_main( hecmesh, fstrsolid, fstrresult, &
          & fstrsolid%SHELL%LAYER(i)%PLUS,  nitem, iitem, ncomp, eitem, jitem, ecomp, i+1, clyr(2*i-1) )
        call fstr_make_result_main( hecmesh, fstrsolid, fstrresult, &
          & fstrsolid%SHELL%LAYER(i)%MINUS, nitem, iitem, ncomp, eitem, jitem, ecomp, i+1, clyr(2*i  ) )
      enddo
      deallocate(clyr)
    endif
 
    ! --- THERMAL STRAIN @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(12) .and. associated(tnstrain) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(12), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'THERMAL_NodalSTRAIN'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = tnstrain(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- CONTACT NORMAL FORCE @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(30) .and. associated(fstrsolid%CONT_NFORCE) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(30), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'CONTACT_NFORCE'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%CONT_NFORCE(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- CONTACT FRICTION FORCE @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(31) .and. associated(fstrsolid%CONT_FRIC) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(31), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'CONTACT_FRICTION'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%CONT_FRIC(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- CONTACT RELATIVE VELOCITY @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(32) .and. associated(fstrsolid%CONT_RELVEL) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(32), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'CONTACT_RELVEL'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%CONT_RELVEL(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- CONTACT STATE @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(33) .and. associated(fstrsolid%CONT_STATE) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(33), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'CONTACT_STATE'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%CONT_STATE(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- CONTACT NORMAL TRACTION @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(36) .and. associated(fstrsolid%CONT_NTRAC) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(36), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'CONTACT_NTRACTION'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%CONT_NTRAC(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- CONTACT FRICTION TRACTION @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(37) .and. associated(fstrsolid%CONT_FTRAC) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(37), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'CONTACT_FTRACTION'
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = fstrsolid%CONT_FTRAC(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- NODE ID @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(38) ) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(38), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'NODE_ID'
      do i = 1, hecmesh%n_node
        fstrresult%node_val_item(nitem*(i-1)+1+iitem) = hecmesh%global_node_ID(i)
      enddo
      iitem = iitem + nn
    endif
 
    ! --- TEMPERATURE @node
    !if( fstrSOLID%output_ctrl(4)%outinfo%on(41) .and. associated(fstrSOLID%CONT_FTRAC) ) then
    !  ncomp = ncomp + 1
    !  nn = n_comp_valtype( fstrSOLID%output_ctrl(4)%outinfo%vtype(41), ndof )
    !  fstrRESULT%nn_dof(ncomp) = nn
    !  fstrRESULT%node_label(ncomp) = 'TEMPERATURE'
    !  do i = 1, hecMESH%n_node
    !    fstrRESULT%node_val_item(nitem*(i-1)+1+iitem) = fstrSOLID%global_node_ID(i)
    !  enddo
    !  iitem = iitem + nn
    !endif
 
 
    ! --- MATERIAL @elem
    if(fstrsolid%output_ctrl(4)%outinfo%on(34)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(34), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'MATERIAL_ID'
      ! hecMESH%section_ID holds the MPC or contact pair index for the pseudo
      ! elements of HEC-MW, not a section ID
      do i = 1, hecmesh%n_elem
        j = hecmesh%section_ID(i)
        if( j >= 1 .and. j <= hecmesh%section%n_sect ) then
          fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = hecmesh%section%sect_mat_ID_item(j)
        else
          fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = 0.d0
        endif
      enddo
      jitem = jitem + nn
    endif
 
    ! --- ELEM ID @elem
    if(fstrsolid%output_ctrl(4)%outinfo%on(39)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(39), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ELEM_ID'
      do i = 1, hecmesh%n_elem
        fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = hecmesh%global_elem_ID(i)
      enddo
      jitem = jitem + nn
    endif
 
    ! --- SECTION ID @elem
    if(fstrsolid%output_ctrl(4)%outinfo%on(40)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(40), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'SECTION_ID'
      do i = 1, hecmesh%n_elem
        fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = hecmesh%section_ID(i)
      enddo
      jitem = jitem + nn
    endif
 
    ! --- ELEMACT flag @element
    if( fstrsolid%output_ctrl(4)%outinfo%on(44) ) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(44), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ELEMACT'
      do i = 1, hecmesh%n_elem
        if( fstrsolid%elements(i)%elemact_flag /= kelact_inactive ) then
          fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = 1.d0
        else
          fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = 0.d0
        end if
      enddo
      jitem = jitem + nn
    endif
 
    if( present(label_suffix) ) then
      do i = 1, fstrresult%nn_component
        fstrresult%node_label(i) = trim(fstrresult%node_label(i))//trim(label_suffix)
      enddo
      do i = 1, fstrresult%ne_component
        fstrresult%elem_label(i) = trim(fstrresult%elem_label(i))//trim(label_suffix)
      enddo
    endif
 
  end subroutine fstr_make_result
 
  !C***
  !C***
  subroutine fstr_merge_result( hecMESH, fstrRESULT, fstrRESULT_add )
    use m_fstr
    use hecmw_util
 
    implicit none
    type (hecmwst_local_mesh) :: hecmesh
    type (hecmwst_result_data) :: fstrresult
    type (hecmwst_result_data) :: fstrresult_add
 
    type (hecmwst_result_data) :: res
    integer(kind=kint) :: i, j, gitem, nitem, eitem, nitem_a, eitem_a
 
    gitem   = sum( fstrresult%ng_dof(1:fstrresult%ng_component) )
    nitem   = sum( fstrresult%nn_dof(1:fstrresult%nn_component) )
    eitem   = sum( fstrresult%ne_dof(1:fstrresult%ne_component) )
    nitem_a = sum( fstrresult_add%nn_dof(1:fstrresult_add%nn_component) )
    eitem_a = sum( fstrresult_add%ne_dof(1:fstrresult_add%ne_component) )
 
    call hecmw_nullify_result_data( res )
    res%ng_component = fstrresult%ng_component
    res%nn_component = fstrresult%nn_component + fstrresult_add%nn_component
    res%ne_component = fstrresult%ne_component + fstrresult_add%ne_component
    allocate( res%ng_dof(res%ng_component), res%global_label(res%ng_component) )
    allocate( res%nn_dof(res%nn_component), res%node_label(res%nn_component) )
    allocate( res%ne_dof(res%ne_component), res%elem_label(res%ne_component) )
    allocate( res%global_val_item(gitem) )
    allocate( res%node_val_item((nitem+nitem_a)*hecmesh%n_node) )
    allocate( res%elem_val_item((eitem+eitem_a)*hecmesh%n_elem) )
 
    res%ng_dof(1:res%ng_component) = fstrresult%ng_dof(1:res%ng_component)
    res%global_label(1:res%ng_component) = fstrresult%global_label(1:res%ng_component)
    res%global_val_item(1:gitem) = fstrresult%global_val_item(1:gitem)
 
    res%nn_dof(1:fstrresult%nn_component) = fstrresult%nn_dof(1:fstrresult%nn_component)
    res%nn_dof(fstrresult%nn_component+1:) = fstrresult_add%nn_dof(1:fstrresult_add%nn_component)
    res%node_label(1:fstrresult%nn_component) = fstrresult%node_label(1:fstrresult%nn_component)
    res%node_label(fstrresult%nn_component+1:) = fstrresult_add%node_label(1:fstrresult_add%nn_component)
    do i = 1, hecmesh%n_node
      do j = 1, nitem
        res%node_val_item((nitem+nitem_a)*(i-1)+j) = fstrresult%node_val_item(nitem*(i-1)+j)
      enddo
      do j = 1, nitem_a
        res%node_val_item((nitem+nitem_a)*(i-1)+nitem+j) = fstrresult_add%node_val_item(nitem_a*(i-1)+j)
      enddo
    enddo
 
    res%ne_dof(1:fstrresult%ne_component) = fstrresult%ne_dof(1:fstrresult%ne_component)
    res%ne_dof(fstrresult%ne_component+1:) = fstrresult_add%ne_dof(1:fstrresult_add%ne_component)
    res%elem_label(1:fstrresult%ne_component) = fstrresult%elem_label(1:fstrresult%ne_component)
    res%elem_label(fstrresult%ne_component+1:) = fstrresult_add%elem_label(1:fstrresult_add%ne_component)
    do i = 1, hecmesh%n_elem
      do j = 1, eitem
        res%elem_val_item((eitem+eitem_a)*(i-1)+j) = fstrresult%elem_val_item(eitem*(i-1)+j)
      enddo
      do j = 1, eitem_a
        res%elem_val_item((eitem+eitem_a)*(i-1)+eitem+j) = fstrresult_add%elem_val_item(eitem_a*(i-1)+j)
      enddo
    enddo
 
    call hecmw_result_free( fstrresult )
    fstrresult = res
  end subroutine fstr_merge_result
 
  subroutine fstr_make_result_main( hecMESH, fstrSOLID, fstrRESULT, RES, nitem, &
      &                              iitem, ncomp, eitem, jitem, ecomp, nlyr, clyr )
    use m_fstr
    use m_out
    use m_static_lib
    use mmaterial
    use hecmw_util
 
    implicit none
    type (hecmwST_local_mesh) :: hecMESH
    type (fstr_solid)         :: fstrSOLID
    type (hecmwST_result_data):: fstrRESULT
    type (fstr_solid_physic_val) :: RES
    integer(kind=kint)        :: istep, flag
    integer(kind=kint)        :: n_lyr, cid
 
    character(len=HECMW_HEADER_LEN) :: header
    character(len=HECMW_NAME_LEN)   :: s, label, nameID, addfname
    character(len=16)                :: clyr
    character(len=12)                :: cnum
    integer(kind=kint) :: i, j, k, ndof, mdof, id, nitem, eitem, nn, mm, ngauss, it
    integer(kind=kint) :: iitem, ncomp, jitem, ecomp, nlyr
    integer(kind=kint) :: isect, cdsys_ID, serr
    real(kind=kreal) :: coords(3,3), coordsys(3, 3)
    ndof = hecmesh%n_dof
 
    ! --- STRAIN @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(3)) then
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(3), ndof )
      ncomp = ncomp + 1
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'NodalSTRAIN'//trim(clyr)
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = res%STRAIN(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- STRESS @node
    if(fstrsolid%output_ctrl(4)%outinfo%on(4)) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(4), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'NodalSTRESS'//trim(clyr)
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = res%STRESS((nn)*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- MISES @node
    if(fstrsolid%output_ctrl(4)%outinfo%on(5)) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(5), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'NodalMISES'//trim(clyr)
      do i = 1, hecmesh%n_node
        fstrresult%node_val_item(nitem*(i-1)+1+iitem) = res%MISES(i)
      enddo
      iitem = iitem + nn
    endif
 
    ! --- Princ STRESS @node
    if(fstrsolid%output_ctrl(4)%outinfo%on(19)) then
      ncomp = ncomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(19), ndof )
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'NodalPrincipalSTRESS'//trim(clyr)
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = res%PSTRESS((nn)*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- Princ STRESS Vector @node
    if(fstrsolid%output_ctrl(4)%outinfo%on(23)) then
      do k=1,3
        write(cnum, '(i0)') k
        ncomp = ncomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(23), ndof )
        fstrresult%nn_dof(ncomp) = nn
        fstrresult%node_label(ncomp) = 'NodalPrincipalSTRESSVector'//trim(cnum)//trim(clyr)
        do i = 1, hecmesh%n_node
          do j = 1, nn
            fstrresult%node_val_item(nitem*(i-1)+j+iitem) = res%PSTRESS_VECT((nn)*(i-1)+j,k)
          enddo
        enddo
        iitem = iitem + nn
      end do
    endif
 
    ! --- Princ STRAIN @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(21)) then
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(21), ndof )
      ncomp = ncomp + 1
      fstrresult%nn_dof(ncomp) = nn
      fstrresult%node_label(ncomp) = 'NodalPrincipalSTRAIN'//trim(clyr)
      do i = 1, hecmesh%n_node
        do j = 1, nn
          fstrresult%node_val_item(nitem*(i-1)+j+iitem) = res%PSTRAIN(nn*(i-1)+j)
        enddo
      enddo
      iitem = iitem + nn
    endif
 
    ! --- Princ STRAIN Vector @node
    if( fstrsolid%output_ctrl(4)%outinfo%on(25)) then
      do k=1,3
        write(cnum, '(i0)') k
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(25), ndof )
        ncomp = ncomp + 1
        fstrresult%nn_dof(ncomp) = nn
        fstrresult%node_label(ncomp) = 'NodalPrincipalSTRAINVector'//trim(cnum)//trim(clyr)
        do i = 1, hecmesh%n_node
          do j = 1, nn
            fstrresult%node_val_item(nitem*(i-1)+j+iitem) = res%PSTRAIN_VECT(nn*(i-1)+j,k)
          enddo
        enddo
        iitem = iitem + nn
      enddo
    endif
 
    ! --- STRAIN @elem
    if( fstrsolid%output_ctrl(4)%outinfo%on(6)) then
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(6), ndof )
      ecomp = ecomp + 1
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ElementalSTRAIN'
      do i = 1, hecmesh%n_elem
        do j = 1, nn
          fstrresult%elem_val_item(eitem*(i-1)+j+jitem) = res%ESTRAIN(nn*(i-1)+j)
        enddo
      enddo
      jitem = jitem + nn
    endif
 
    ! --- STRESS @elem
    if(fstrsolid%output_ctrl(4)%outinfo%on(7)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(7), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ElementalSTRESS'
      do i = 1, hecmesh%n_elem
        do j = 1, nn
          fstrresult%elem_val_item(eitem*(i-1)+j+jitem) = res%ESTRESS((nn)*(i-1)+j)
        enddo
      enddo
      jitem = jitem + nn
    endif
 
    ! --- MISES @elem
    if(fstrsolid%output_ctrl(4)%outinfo%on(8)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(8), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ElementalMISES'
      do i = 1, hecmesh%n_elem
        fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = res%EMISES(i)
      enddo
      jitem = jitem + nn
    endif
 
    ! --- Principal_STRESS @element
    if(fstrsolid%output_ctrl(4)%outinfo%on(20)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(20), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ElementalPrincipalSTRESS'
      do i = 1, hecmesh%n_elem
        do j = 1, nn
          fstrresult%elem_val_item(eitem*(i-1)+j+jitem) = res%EPSTRESS((nn)*(i-1)+j)
        enddo
      enddo
      jitem = jitem + nn
    endif
 
    ! --- Principal_STRAIN @element
    if(fstrsolid%output_ctrl(4)%outinfo%on(22)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(22), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ElementalPrincipalSTRAIN'
      do i = 1, hecmesh%n_elem
        do j = 1, nn
          fstrresult%elem_val_item(eitem*(i-1)+j+jitem) = res%EPSTRAIN((nn)*(i-1)+j)
        enddo
      enddo
      jitem = jitem + nn
    endif
 
    ! --- ELEM PRINC STRESS VECTOR
    if(fstrsolid%output_ctrl(4)%outinfo%on(24)) then
      do k = 1, 3
        write(cnum,'(i0)')k
        ecomp = ecomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(24), ndof )
        fstrresult%ne_dof(ecomp) = nn
        fstrresult%elem_label(ecomp) = 'ElementalPrincipalSTRESSVector'//trim(cnum)
        do i = 1, hecmesh%n_elem
          do j = 1, nn
            fstrresult%elem_val_item(eitem*(i-1)+j+jitem) = res%EPSTRESS_VECT((nn)*(i-1)+j,k)
          enddo
        enddo
        jitem = jitem + nn
      enddo
    endif
 
    ! --- ELEM PRINC STRAIN VECTOR
    if(fstrsolid%output_ctrl(4)%outinfo%on(26)) then
      do k = 1, 3
        write(cnum,'(i0)')k
        ecomp = ecomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(26), ndof )
        fstrresult%ne_dof(ecomp) = nn
        fstrresult%elem_label(ecomp) = 'ElementalPrincipalSTRAINVector'//trim(cnum)
        do i = 1, hecmesh%n_elem
          do j = 1, nn
            fstrresult%elem_val_item(eitem*(i-1)+j+jitem) = res%EPSTRAIN_VECT((nn)*(i-1)+j,k)
          enddo
        enddo
        jitem = jitem + nn
      enddo
    endif
 
    ! --- PLSTRAIN @elem
    if(fstrsolid%output_ctrl(4)%outinfo%on(43)) then
      ecomp = ecomp + 1
      nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(43), ndof )
      fstrresult%ne_dof(ecomp) = nn
      fstrresult%elem_label(ecomp) = 'ElementalPLSTRAIN'
      do i = 1, hecmesh%n_elem
        res%EPLSTRAIN(i) = 0.d0
        do j = 1, size(fstrsolid%elements(i)%gausses) 
          res%EPLSTRAIN(i) = res%EPLSTRAIN(i) + fstrsolid%elements(i)%gausses(j)%plstrain
        enddo
        res%EPLSTRAIN(i) = res%EPLSTRAIN(i) / size(fstrsolid%elements(i)%gausses)
        fstrresult%elem_val_item(eitem*(i-1)+1+jitem) = res%EPLSTRAIN(i)
      enddo
      jitem = jitem + nn
    endif
 
    ! --- ORIENTATION @elem
    if(fstrsolid%output_ctrl(4)%outinfo%on(42)) then
      do i = 1, 3
        write(cnum,'(i0)')i
        ecomp = ecomp + 1
        nn = n_comp_valtype( fstrsolid%output_ctrl(4)%outinfo%vtype(42), ndof )
        fstrresult%ne_dof(ecomp) = nn
        fstrresult%elem_label(ecomp) = 'ORIENTATION'//trim(cnum)
      enddo
      ! hecMESH%section_ID holds the MPC or contact pair index for the pseudo
      ! elements of HEC-MW, not a section ID
      do i = 1, hecmesh%n_elem
        isect = hecmesh%section_ID(i)
        if( isect >= 1 .and. isect <= hecmesh%section%n_sect ) then
          cdsys_id = hecmesh%section%sect_orien_ID(isect)
        else
          cdsys_id = -1
        endif
        if(cdsys_id == -1) then
          coordsys(:,:) = 0.0d0
        else
          call get_coordsys(cdsys_id, hecmesh, fstrsolid, coords, i)
          call set_localcoordsys( coords, g_localcoordsys(cdsys_id), coordsys(:, :), serr )
        endif
        do j = 1, 3
          do k = 1, 3
            fstrresult%elem_val_item(eitem*(i-1)+jitem+nn*(j-1)+k) = coordsys(j,k)
          enddo
        enddo
      enddo
      jitem = jitem + 3*nn
    endif
 
  end subroutine fstr_make_result_main
 
  subroutine fstr_reorder_node_shell(fstrSOLID, hecMESH, unode)
    use m_fstr
    use m_out
    use m_static_lib
 
    implicit none
    type (fstr_solid)         :: fstrsolid
    type (hecmwst_local_mesh) :: hecmesh
    integer(kind=kint) :: i, j, k, itype, is, ie, ic_type, js, icel
    integer(kind=kint) :: mm, n1, n2
    real(kind=kreal), allocatable   :: unode(:)
 
    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)
      if(ic_type == 781)then
        do icel = is, ie
          js = hecmesh%elem_node_index(icel-1)
          do j = 1, 4
            n1 = hecmesh%elem_node_item(js+j  )
            n2 = hecmesh%elem_node_item(js+j+4)
            unode(3*n2-2) = unode(3*n1-2)
            unode(3*n2-1) = unode(3*n1-1)
            unode(3*n2  ) = unode(3*n1  )
          enddo
        enddo
      elseif(ic_type == 761)then
        do icel = is, ie
          js = hecmesh%elem_node_index(icel-1)
          do j = 1, 3
            n1 = hecmesh%elem_node_item(js+j  )
            n2 = hecmesh%elem_node_item(js+j+3)
            unode(3*n2-2) = unode(3*n1-2)
            unode(3*n2-1) = unode(3*n1-1)
            unode(3*n2  ) = unode(3*n1  )
          enddo
        enddo
      endif
    enddo
 
  end subroutine fstr_reorder_node_shell
 
  subroutine fstr_reorder_rot_shell(fstrSOLID, hecMESH, unode)
    use m_fstr
    use m_out
    use m_static_lib
 
    implicit none
    type (fstr_solid)         :: fstrsolid
    type (hecmwst_local_mesh) :: hecmesh
    integer(kind=kint) :: i, j, k, itype, is, ie, ic_type, js, icel
    integer(kind=kint) :: mm, n1, n2
    real(kind=kreal), allocatable   :: unode(:)
 
    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)
      if(ic_type == 781)then
        do icel = is, ie
          js = hecmesh%elem_node_index(icel-1)
          do j = 1, 4
            n1 = hecmesh%elem_node_item(js+j)
            n2 = hecmesh%elem_node_item(js+j+4)
            unode(3*n1-2) = fstrsolid%unode(3*n2-2)
            unode(3*n1-1) = fstrsolid%unode(3*n2-1)
            unode(3*n1  ) = fstrsolid%unode(3*n2  )
            unode(3*n2-2) = fstrsolid%unode(3*n2-2)
            unode(3*n2-1) = fstrsolid%unode(3*n2-1)
            unode(3*n2  ) = fstrsolid%unode(3*n2  )
          enddo
        enddo
      elseif(ic_type == 761)then
        do icel = is, ie
          js = hecmesh%elem_node_index(icel-1)
          do j = 1, 3
            n1 = hecmesh%elem_node_item(js+j)
            n2 = hecmesh%elem_node_item(js+j+3)
 
            unode(3*n1-2) = fstrsolid%unode(3*n2-2)
            unode(3*n1-1) = fstrsolid%unode(3*n2-1)
            unode(3*n1  ) = fstrsolid%unode(3*n2  )
            unode(3*n2-2) = fstrsolid%unode(3*n2-2)
            unode(3*n2-1) = fstrsolid%unode(3*n2-1)
            unode(3*n2  ) = fstrsolid%unode(3*n2  )
          enddo
        enddo
      endif
    enddo
 
  end subroutine fstr_reorder_rot_shell
 
  subroutine fstr_reorder_node_beam(fstrSOLID, hecMESH, unode)
    use m_fstr
    use m_out
    use m_static_lib
 
    implicit none
    type (fstr_solid)         :: fstrsolid
    type (hecmwst_local_mesh) :: hecmesh
    integer(kind=kint) :: i, j, k, itype, is, ie, ic_type, js, icel
    integer(kind=kint) :: mm, a, b
    real(kind=kreal), allocatable   :: unode(:)
 
    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)
      if(ic_type == 641)then
        do icel = is, ie
          js = hecmesh%elem_node_index(icel-1)
          do j = 1, 2
            a = hecmesh%elem_node_item(js+j)
            b = hecmesh%elem_node_item(js+j+2)
            unode(3*b-2) = unode(3*a-2)
            unode(3*b-1) = unode(3*a-1)
            unode(3*b  ) = unode(3*a  )
          enddo
        enddo
      endif
    enddo
 
  end subroutine fstr_reorder_node_beam
 
  subroutine setup_contact_output_variables( hecMESH, fstrSOLID, phase, dtime, expflag )
    use m_fstr
    use hecmw_util
    use mcontact
    implicit none
    type(hecmwst_local_mesh), intent(in)  :: hecmesh
    type (fstr_solid), intent(inout)      :: fstrsolid
    integer(kind=kint), intent(in)        :: phase
    real(kind=kreal), optional, intent(in) :: dtime   
    logical, optional, intent(in) :: expflag
 
    integer(kind=kint), parameter :: nval = 10
    logical, save :: updated(nval) = .false.
    integer(kind=kint) :: ndof, i
    real(kind=kreal) :: area, dt_use
 
    ndof = hecmesh%n_dof
 
    if( phase == -1 ) then
      updated(1:nval) = .false.
      return
    else
      if( phase /= 3 .and. phase /= 4 ) return !irregular case
    end if
 
    ! Update contact state vectors if dtime is provided
    if( present(dtime) .and. dtime > 0.d0 ) then
      dt_use = dtime
    else
      dt_use = 1.d0  ! default value if dtime not available
    endif
    if( .not. updated(3) .and. .not. updated(4) ) then
      call fstr_update_contact_state_vectors( fstrsolid, dt_use )
    endif
 
    ! --- CONTACT NORMAL FORCE @node
    if( fstrsolid%output_ctrl(phase)%outinfo%on(30) .and. associated(fstrsolid%CONT_NFORCE) ) then
      if( .not. updated(1) ) then
        ! Explicit dynamic has no Lagrange matrix, so the solver cannot fill CONT_NFORCE
        ! via fstr_calc_contact_output_force. Build it here from the stored contact
        ! multipliers before the parallel aggregation. Implicit/static keep the value
        ! already computed inside the solver (flag absent), so this is a no-op for them.
        if( present(expflag) .and. expflag ) then
          call fstr_calc_contact_output_force_exp( hecmesh, fstrsolid )
        end if
        if( paracontactflag ) then
          call fstr_setup_parancon_contactvalue(hecmesh,ndof,fstrsolid%CONT_NFORCE,1)
        end if
      end if
      updated(1) = .true.
    endif
 
    ! --- CONTACT FRICTION FORCE @node
    if( fstrsolid%output_ctrl(phase)%outinfo%on(31) .and. associated(fstrsolid%CONT_FRIC) ) then
      if( paracontactflag .and. .not. updated(2)) then
        call fstr_setup_parancon_contactvalue(hecmesh,ndof,fstrsolid%CONT_FRIC,1)
      end if
      updated(2) = .true.
    endif
 
    ! --- CONTACT RELATIVE VELOCITY @node
    if( fstrsolid%output_ctrl(phase)%outinfo%on(32) .and. associated(fstrsolid%CONT_RELVEL) ) then
      if( paracontactflag .and. .not. updated(3)) then
        call fstr_setup_parancon_contactvalue(hecmesh,ndof,fstrsolid%CONT_RELVEL,1)
      end if
      updated(3) = .true.
    endif
 
    ! --- CONTACT STATE @node
    if( fstrsolid%output_ctrl(phase)%outinfo%on(33) .and. associated(fstrsolid%CONT_STATE) ) then
      if( paracontactflag .and. .not. updated(4)) then
        call fstr_setup_parancon_contactvalue(hecmesh,1,fstrsolid%CONT_STATE,2)
      end if
      updated(4) = .true.
    endif
 
    ! --- CONTACT AREA for CONTACT TRACTION
    if( fstrsolid%output_ctrl(phase)%outinfo%on(36) .or. fstrsolid%output_ctrl(phase)%outinfo%on(37) ) then
      if( .not. updated(5)) call calc_contact_area( hecmesh, fstrsolid, 0 )
      ! fstr_setup_parancon_contactvalue is not necessary because
      ! contact area is calculated from original surface group
    end if
 
    ! --- CONTACT NORMAL TRACTION @node
    if( fstrsolid%output_ctrl(phase)%outinfo%on(36) .and. associated(fstrsolid%CONT_NTRAC) ) then
      if( paracontactflag .and. .not. updated(6)) then
        if( .not. updated(1)) call fstr_setup_parancon_contactvalue(hecmesh,ndof,fstrsolid%CONT_NFORCE,1)
      end if
      fstrsolid%CONT_NTRAC(:) = 0.d0
      do i=1,hecmesh%nn_internal
        area = fstrsolid%CONT_AREA(i)
        if( area < 1.d-16 ) cycle
        fstrsolid%CONT_NTRAC(3*i-2:3*i) = fstrsolid%CONT_NFORCE(3*i-2:3*i)/area
      end do
      updated(6) = .true.
    endif
 
    ! --- CONTACT FRICTION TRACTION @node
    if( fstrsolid%output_ctrl(phase)%outinfo%on(37) .and. associated(fstrsolid%CONT_FTRAC) ) then
      if( paracontactflag .and. .not. updated(7)) then
        if( .not. updated(1)) call fstr_setup_parancon_contactvalue(hecmesh,ndof,fstrsolid%CONT_FRIC,1)
      end if
      fstrsolid%CONT_FTRAC(:) = 0.d0
      do i=1,hecmesh%nn_internal
        area = fstrsolid%CONT_AREA(i)
        if( area < 1.d-16 ) cycle
        fstrsolid%CONT_FTRAC(3*i-2:3*i) = fstrsolid%CONT_FRIC(3*i-2:3*i)/area
      end do
      updated(7) = .true.
    endif
 
  end subroutine
 
end module m_make_result