Dear Friends,
I have a simple UEL subroutine with a 2D truss element. The element has two nodes: 11 (down) and 12 (up). The nodes are connected from bottom to top. Boundary conditions: node 11 has pinned support; node 12 has a roller (x-direction is blocked). The truss element is loaded by applying displacement U in the positive Y-direction at node 12.
Below, you can find my UEL subroutine and a sketch of the problem:
The issue is, that the sign of calculated in UEL reaction forces do not coincide with results from .dat-File.
If I understand correctly, the RHS for this problem is equivalent to the internal forces of a truss element in the global coordinate system. I may be wrong, but I think that the internal forces of a truss element in the global coordinate system are equal to the reaction forces of the truss element.
To calculate the RHS, I follow this procedure:
- Define the global displacement vector U_glob. This vector is equal to the displacement vector U [U(1), U(2), U(3), U(4)] which is calculated by Abaqus.
- Define the local displacement vector U_loc. U_loc = T*U_glob. T is a transformation matrix.
- Define the local force vector localF. localF = K*U_loc. K is an element stiffness matrix in the local coordinate system.
- Define the global force vector globalF. globalF = tranT * localF. tranT is the transpose of matrix T.
- RHS = globalF
As you can see in the output of the Abaqus commander, the reaction force at node 11 (RF11) is negative, and the reaction force at the top node 12 (N) is positive. Considering the direction of the applied loading U, the sign of the result is correct. However, in the dat-file, the reaction force at node 11 (RF2) is positive and at node 12 it's negative. This is an incorrect result.
I also have a reference element T2D2 (nodes 1-2) in Abaqus CAE to verify the results. As you can see in output from dat-File, the sign is correct.
Do you have any suggestions to resolve my question?
Many thanks in advance!
============================================================== C The Fortran code below represents a User Defined Truss Element SUBROUTINE UEL(RHS,AMATRX,SVARS,ENERGY,NDOFEL,NRHS,NSVARS, 1 PROPS,NPROPS,COORDS,MCRD,NNODE,U,DU,V,A,JTYPE,TIME, 2 DTIME,KSTEP,KINC,JELEM,PARAMS,NDLOAD,JDLTYP,ADLMAG, 3 PREDEF,NPREDF,LFLAGS,MLVARX,DDLMAG,MDLOAD,PNEWDT, 4 JPROPS,NJPROP,PERIOD) INCLUDE 'ABA_PARAM.INC' DIMENSION RHS(NDOFEL,NRHS),AMATRX(NDOFEL,NDOFEL), 1 SVARS(NSVARS),ENERGY(8),PROPS(*),COORDS(MCRD,NNODE), 2 U(NDOFEL),DU(MLVARX,*),V(NDOFEL),A(NDOFEL),TIME(2), 3 PARAMS(3),JDLTYP(MDLOAD,*),ADLMAG(MDLOAD,*), 4 DDLMAG(MDLOAD,*),PREDEF(2,NPREDF,NNODE),LFLAGS(*), 5 JPROPS(*) C This is a UEL Subroutine for a Truss Element integer i, j, m integer, parameter :: ikind=selected_real_kind(p=8) real(kind=ikind):: COSINEBETA, SINEBETA, L, F real(kind=ikind):: K(2,2), T(2,4), tranT(4,2) real(kind=ikind):: midK(4,2), globalK(4,4) real(kind=ikind):: localF(2,1), globalF(4,1), U_loc(2,1), U_glob(4,1) E = PROPS(1) AREA = PROPS(2) x1 = COORDS(1,1) write(*,*) 'x1 =', x1 y1 = COORDS(2,1) write(*,*) 'y1 =', y1 x2 = COORDS(1,2) write(*,*) 'x2 =', x2 y2 = COORDS(2,2) write(*,*) 'y2 =', y2 write(*,*) 'COORDS =', COORDS C Calculating the length of the truss element L = sqrt((x2-x1)**2 + (y2-y1)**2) COSINEBETA = (x2-x1)/L SINEBETA = (y2-y1)/L C Zeroing the stiffness matrix DO K1 = 1, NDOFEL DO KRHS = 1, NRHS RHS(K1,KRHS) = ZERO END DO DO K2 = 1, NDOFEL AMATRX(K2,K1) = ZERO END DO END DO C Defining a local stiffness matrix K(1,1) = E*AREA/L K(1,2) = -E*AREA/L K(2,1) = -E*AREA/L K(2,2) = E*AREA/L C Printing the local stiffness matrix write(*,*) 'local stiffness matri: K =' do i = 1, 2 print *,(K(i,j), j = 1,2) end do write(*,*)'' C Defining a transformation matrix T(1,1) = COSINEBETA T(1,2) = SINEBETA T(1,3) = 0 T(1,4) = 0 T(2,1) = 0 T(2,2) = 0 T(2,3) = COSINEBETA T(2,4) = SINEBETA C Printing a transformation matrix write(*,*) 'transformation matrix: T =' do i = 1, 2 print *,(T(i,j), j = 1,4) end do write(*,*)'' C Calculation of the transpose of a transformation matrix do i=1, 2 do j=1, 4 tranT(j,i) = T(i,j) end do end do write(*,*) 'transpose of a transformation matrix T: tranT =' do i = 1, 4 print *,(tranT(i,j), j = 1,2) end do write(*,*)'' C midK is the multiplication of tranT and K, so tranT*K do i=1, 4 do j=1, 2 midK(i, j) = 0.0 do m=1, 2 midK(i,j) = midK(i,j) + tranT(i, m)*K(m, j) end do end do end do write(*,*) 'multiplication of tranT and K, so tranT*K =' do i = 1, 4 print *,(midK(i,j), j = 1,2) end do write(*,*)'' C globalK is the global stiffness matrix of a truss element do i=1, 4 do j=1, 4 globalK(i,j) = 0.0 do m = 1, 2 globalK(i,j) = globalK(i,j) + midK(i, m) * T(m, j) end do end do end do write(*,*) 'global stiffness matrix of a truss element, globalK =' do i = 1, 4 print *,(globalK(i,j), j = 1,4) end do write(*,*)'' C Normal implicit time incrementation procedure IF (LFLAGS(3) .EQ. 1) THEN C Static automatic and direct incrementation if(LFLAGS(1) .EQ. 1 .OR. LFLAGS(1) .EQ. 2) THEN AMATRX = globalK write(*,*) 'AMATRX =' do i = 1, 4 print *,(AMATRX(i,j), j = 1,4) end do write(*,*)'' C Defining the global displacement vector U_glob(1,1) = U(1) U_glob(2,1) = U(2) U_glob(3,1) = U(3) U_glob(4,1) = U(4) write(*,*) 'U_glob=' ! Printing the resulting U_glob matrix DO I = 1, 4 WRITE(*,*) U_glob(I, 1) END DO C Defining the local displacement vector do i = 1, 2 U_loc(i, 1) = 0.0 do j = 1, 4 U_loc(i, 1) = U_loc(i, 1) + T(i, j) * U_glob(j, 1) end do end do ! Printing the resulting U_loc matrix write(*,*) 'U_loc=' do i = 1, 2 write(*,*) U_loc(i, 1) end do C Defining the local and global force vector localF(1,1) = K(1,1)*U_loc(1,1) + K(1,2)*U_loc(2,1) localF(2,1) = K(2,1)*U_loc(1,1) + K(2,2)*U_loc(2,2) print *, 'localF= ', localF do i=1, 4 globalF(i, 1) = 0.0 do j=1, 2 globalF(i,1) = globalF(i,1) + tranT(i,j)*localF(j, 1) end do end do C Identifying the right hand side force vector RHS = globalF print *, 'RHS=globalF=Reaction forces=', RHS end if end if write(*,*)'' write(*,*)'U(2)=', U(2) write(*,*)'U(4)=', U(4) write(*,*)'END OF ONE ITERATION' write(*,*)'********************' write(*,*)'' RETURN END ======================================================== Abaqus