cdefo3.F

Go to the documentation of this file.

Copyright>        OpenRadioss
Copyright>        Copyright (C) 1986-2025 Altair Engineering Inc.
Copyright>
Copyright>        This program is free software: you can redistribute it and/or modify
Copyright>        it under the terms of the GNU Affero General Public License as published by
Copyright>        the Free Software Foundation, either version 3 of the License, or
Copyright>        (at your option) any later version.
Copyright>
Copyright>        This program is distributed in the hope that it will be useful,
Copyright>        but WITHOUT ANY WARRANTY; without even the implied warranty of
Copyright>        MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
Copyright>        GNU Affero General Public License for more details.
Copyright>
Copyright>        You should have received a copy of the GNU Affero General Public License
Copyright>        along with this program.  If not, see <https://www.gnu.org/licenses/>.
Copyright>
Copyright>
Copyright>        Commercial Alternative: Altair Radioss Software
Copyright>
Copyright>        As an alternative to this open-source version, Altair also offers Altair Radioss
Copyright>        software under a commercial license.  Contact Altair to discuss further if the
Copyright>        commercial version may interest you: https://www.altair.com/radioss/.
!||====================================================================
!||    cdefo3       ../engine/source/elements/shell/coque/cdefo3.F
!||--- called by ------------------------------------------------------
!||    cforc3       ../engine/source/elements/shell/coque/cforc3.F
!||    cforc3_crk   ../engine/source/elements/xfem/cforc3_crk.f
!||====================================================================
      SUBROUTINE cdefo3(JFT,JLT,VL1,VL2,VL3,VL4,
     1            DT1C,PX1 ,PX2,PY1,PY2,AREA,
     2            EXX ,EYY,EXY,EXZ,EYZ ,X2  ,X3 ,X4 ,Y2 ,Y3  ,
     3            Y4  ,Z2 ,VX1,VX2,VX3 ,VX4 ,VY1,VY2,VY3,VY4 ,
     4            VZ1 ,VZ2,VZ3,VZ4,E1X ,E1Y ,E1Z,E2X,E2Y,E2Z ,
     5            E3X ,E3Y,E3Z,IHBE)
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   G l o b a l   P a r a m e t e r s
C-----------------------------------------------
#include      "mvsiz_p.inc"
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "impl1_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER JFT, JLT,IHBE
C     REAL
      my_real
     .   VL1(MVSIZ,3),VL2(MVSIZ,3),VL3(MVSIZ,3),VL4(MVSIZ,3),DT1C(*),
     .   PX1(MVSIZ), PX2(MVSIZ), PY1(MVSIZ), PY2(MVSIZ),
     .   exx(mvsiz), area(mvsiz),
     .   eyy(mvsiz), exy(mvsiz), exz(mvsiz), eyz(mvsiz),
     .   x2(mvsiz), x3(mvsiz), x4(mvsiz),
     .   y2(mvsiz), y3(mvsiz), y4(mvsiz), z2(mvsiz),
     .   vx1(mvsiz),  vx2(mvsiz),  vx3(mvsiz),  vx4(mvsiz),
     .   vy1(mvsiz), vy2(mvsiz), vy3(mvsiz), vy4(mvsiz),
     .   vz1(mvsiz), vz2(mvsiz), vz3(mvsiz), vz4(mvsiz),
     .   vx13(mvsiz), vx24(mvsiz), vy13(mvsiz), vy24(mvsiz),
     .   vz13(mvsiz), vz24(mvsiz), 
     .   e1x(mvsiz), e1y(mvsiz), e1z(mvsiz), e2x(mvsiz),
     .   e2y(mvsiz), e2z(mvsiz), e3x(mvsiz), e3y(mvsiz), e3z(mvsiz)
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I, J
      my_real
     .     ZZ2(MVSIZ), EXZZ2(MVSIZ), EYZZ2(MVSIZ), EXZ2PY2(MVSIZ), 
     .     EXZ2PY1(MVSIZ), EYZ2PX2(MVSIZ), EYZ2PX1(MVSIZ),ZZZ(MVSIZ), 
     .     EXZ2(MVSIZ),EYZ2(MVSIZ),
     .     dt1v4, tmp1a(mvsiz), tmp2a(mvsiz), tmp3a(mvsiz),
     .     tmp1b(mvsiz), tmp2b(mvsiz), tmp3b(mvsiz),gzx(mvsiz),
     .     gzy(mvsiz)
C=======================================================================
C
      DO i=jft,jlt
        vx1(i)=e1x(i)*vl1(i,1)+e1y(i)*vl1(i,2)+e1z(i)*vl1(i,3)
        vx2(i)=e1x(i)*vl2(i,1)+e1y(i)*vl2(i,2)+e1z(i)*vl2(i,3)
        vx3(i)=e1x(i)*vl3(i,1)+e1y(i)*vl3(i,2)+e1z(i)*vl3(i,3)
        vx4(i)=e1x(i)*vl4(i,1)+e1y(i)*vl4(i,2)+e1z(i)*vl4(i,3)
C
        vy4(i)=e2x(i)*vl4(i,1)+e2y(i)*vl4(i,2)+e2z(i)*vl4(i,3)
        vy3(i)=e2x(i)*vl3(i,1)+e2y(i)*vl3(i,2)+e2z(i)*vl3(i,3)
        vy2(i)=e2x(i)*vl2(i,1)+e2y(i)*vl2(i,2)+e2z(i)*vl2(i,3)
        vy1(i)=e2x(i)*vl1(i,1)+e2y(i)*vl1(i,2)+e2z(i)*vl1(i,3)
C
        vz1(i)=e3x(i)*vl1(i,1)+e3y(i)*vl1(i,2)+e3z(i)*vl1(i,3)
        vz2(i)=e3x(i)*vl2(i,1)+e3y(i)*vl2(i,2)+e3z(i)*vl2(i,3)
        vz3(i)=e3x(i)*vl3(i,1)+e3y(i)*vl3(i,2)+e3z(i)*vl3(i,3)
        vz4(i)=e3x(i)*vl4(i,1)+e3y(i)*vl4(i,2)+e3z(i)*vl4(i,3)
C
        vz13(i)=vz1(i)-vz3(i)
        vz24(i)=vz2(i)-vz4(i)
C
        eyz(i)=py1(i)*vz13(i)+py2(i)*vz24(i)
        exz(i)=px1(i)*vz13(i)+px2(i)*vz24(i)
C
      ENDDO
C
      IF(ihbe<=1)THEN
        DO i=jft,jlt
C z2 mis a zero pour non equilibre des foreces et moments d'hourglass
          z2(i) = zero
          dt1v4 = fourth*dt1c(i)
          IF(impl_s>0) dt1v4=zero
          tmp2a(i)=py2(i)+py1(i)
          tmp3a(i)=sign(max(abs(tmp2a(i)),em20),tmp2a(i))
          tmp1a(i)=dt1v4*(vz13(i)-vz24(i))*(vz13(i)-vz24(i))/tmp3a(i)
          vx13(i)=vx1(i)-vx3(i)
          vx24(i)=vx2(i)-vx4(i)
          vx13(i)=vx13(i)-tmp1a(i)
          vx24(i)=vx24(i)+tmp1a(i)
C
          exx(i)=px1(i)*vx13(i)+px2(i)*vx24(i)
          exy(i)=py1(i)*vx13(i)+py2(i)*vx24(i)
C
          tmp1b(i)=px2(i)-px1(i)
          tmp3b(i)=sign(max(abs(tmp1b(i)),em20),tmp1b(i))
          tmp2b(i)=dt1v4*(vz13(i)+vz24(i))*(vz13(i)+vz24(i))/tmp3b(i)
          vy13(i)=vy1(i)-vy3(i)
          vy24(i)=vy2(i)-vy4(i)
          vy13(i)=vy13(i)+tmp2b(i)
          vy24(i)=vy24(i)+tmp2b(i)
C
          exy(i)=exy(i)+px1(i)*vy13(i)+px2(i)*vy24(i)
          eyy(i)=py1(i)*vy13(i)+py2(i)*vy24(i)
        ENDDO
C
      ELSEIF(ihbe==2.OR.ihbe==3)THEN
C
        DO i=jft,jlt
          dt1v4 = half*dt1c(i)
          IF(impl_s>0) dt1v4=zero
C
          gzx(i)  = exz(i)/area(i)
          exzz2(i)= gzx(i)*z2(i)
          exz2(i) = gzx(i)*gzx(i)*dt1v4
          vx3(i)=vx3(i)         -exz2(i)*x3(i)-vx1(i)
          vx2(i)=vx2(i)+exzz2(i)-exz2(i)*x2(i)-vx1(i)
          vx4(i)=vx4(i)+exzz2(i)-exz2(i)*x4(i)-vx1(i)
          vx1(i)=zero
C
          gzy(i)  = eyz(i)/area(i)
          eyzz2(i)= gzy(i)*z2(i)
          eyz2(i) = gzy(i)*gzy(i)*dt1v4
          vy3(i)=vy3(i)         -eyz2(i)*y3(i)-vy1(i)
          vy2(i)=vy2(i)+eyzz2(i)-eyz2(i)*y2(i)-vy1(i)
          vy4(i)=vy4(i)+eyzz2(i)-eyz2(i)*y4(i)-vy1(i)
          vy1(i)=zero
C
          zzz(i)= (exz2(i)+eyz2(i))*z2(i)
          vz3(i)=vz3(i)-gzy(i)*y3(i)-gzx(i)*x3(i)       -vz1(i)
          vz2(i)=vz2(i)-gzy(i)*y2(i)-gzx(i)*x2(i)-zzz(i)-vz1(i)
          vz4(i)=vz4(i)-gzy(i)*y4(i)-gzx(i)*x4(i)-zzz(i)-vz1(i)
          vz1(i)=zero
C
          vx13(i)=-vx3(i)
          vx24(i)=vx2(i)-vx4(i)
C
          exx(i)=px1(i)*vx13(i)+px2(i)*vx24(i)
          exy(i)=py1(i)*vx13(i)+py2(i)*vx24(i)
C
          vy13(i)=-vy3(i)
          vy24(i)=vy2(i)-vy4(i)
C
          exy(i)=exy(i)+px1(i)*vy13(i)+px2(i)*vy24(i)
          eyy(i)=py1(i)*vy13(i)+py2(i)*vy24(i)
        ENDDO
C
      ELSEIF(ihbe==4)THEN
C    old
        DO i=jft,jlt
          dt1v4 = half*dt1c(i)
          IF(impl_s>0) dt1v4=zero
c wraping z1=z3=0 z4=z2 => modif des vitesses pour l'hourglass membrane
          zz2(i)=half*z2(i)
C
          gzx(i)     = exz(i)/area(i)
          exzz2(i)   = gzx(i)*zz2(i)
          exz2(i)    = gzx(i)*gzx(i)*dt1v4
          exz2py2(i) = exz2(i)*py2(i)
          exz2py1(i) = exz2(i)*py1(i)
          vx1(i)=vx1(i)-exzz2(i)-exz2py2(i)
          vx3(i)=vx3(i)-exzz2(i)+exz2py2(i)
          vx2(i)=vx2(i)+exzz2(i)+exz2py1(i)
          vx4(i)=vx4(i)+exzz2(i)-exz2py1(i)
C
          gzy(i)     = eyz(i)/area(i)
          eyzz2(i)   = gzy(i)*zz2(i)
          eyz2(i)    = gzy(i)*gzy(i)*dt1v4
          eyz2px2(i) = eyz2(i)*px2(i)
          eyz2px1(i) = eyz2(i)*px1(i)
          vy1(i)=vy1(i)-eyzz2(i)+eyz2px2(i)
          vy3(i)=vy3(i)-eyzz2(i)-eyz2px2(i)
          vy2(i)=vy2(i)+eyzz2(i)-eyz2px1(i)
          vy4(i)=vy4(i)+eyzz2(i)+eyz2px1(i)
C
          vx13(i)=vx1(i)-vx3(i)
          vx24(i)=vx2(i)-vx4(i)
C
          exx(i)=px1(i)*vx13(i)+px2(i)*vx24(i)
          exy(i)=py1(i)*vx13(i)+py2(i)*vx24(i)
C
          vy13(i)=vy1(i)-vy3(i)
          vy24(i)=vy2(i)-vy4(i)
C
          exy(i)=exy(i)+px1(i)*vy13(i)+px2(i)*vy24(i)
          eyy(i)=py1(i)*vy13(i)+py2(i)*vy24(i)
        ENDDO
      ENDIF
c-----------
      RETURN
      END
c
!||====================================================================
!||    cdefot3   ../engine/source/elements/shell/coque/cdefo3.F
!||--- called by ------------------------------------------------------
!||    cforc3    ../engine/source/elements/shell/coque/cforc3.F
!||====================================================================
      SUBROUTINE cdefot3(JFT ,JLT  ,PX1  ,PX2  ,PY1   ,
     2                   PY2 ,VX13 ,VX24  ,VY13  ,VY24   ,
     3                   F_DEF )
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   G l o b a l   P a r a m e t e r s
C-----------------------------------------------
#include      "mvsiz_p.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER JFT, JLT
      my_real
     .   PX1(MVSIZ), PX2(MVSIZ), PY1(MVSIZ), PY2(MVSIZ),
     .   VX13(MVSIZ), VX24(MVSIZ), VY13(MVSIZ), VY24(MVSIZ),
     .   F_DEF(MVSIZ,8)
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I, J
      my_real
     .     S
C=======================================================================
        DO I=jft,jlt
          f_def(i,1)=px1(i)*vx13(i)+px2(i)*vx24(i)
          f_def(i,2)=py1(i)*vy13(i)+py2(i)*vy24(i)
          f_def(i,3)=py1(i)*vx13(i)+py2(i)*vx24(i)
          f_def(i,4)=px1(i)*vy13(i)+px2(i)*vy24(i)
        ENDDO
c-----------
      RETURN
      END