m1ismstr10_pon.F File Reference

#include "implicit_f.inc"
#include "mvsiz_p.inc"
#include "vectorize.inc"

Go to the source code of this file.

Functions/Subroutines
subroutine	m1ismstr10_pon (nel, es1, es2, es3, es4, es5, es6, epsth, c1, g, sig)

Function/Subroutine Documentation

m1ismstr10_pon()

subroutine m1ismstr10_pon	(	integer	nel,
			es1,
			es2,
			es3,
			es4,
			es5,
			es6,
			epsth,
			c1,
			g,
			sig )

Definition at line 30 of file m1ismstr10_pon.F.

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 
C-----------------------------------------------
C----------------------------------------------------------------
C  I N P U T   A R G U M E N T S
C----------------------------------------------------------------
      INTEGER       NEL     
      my_real
     .   es1(*), es2(*), es3(*), 
     .   es4(*), es5(*), es6(*),
     .   sig(mvsiz,6), epsth(*),c1(*),g(*)
C----------------------------------------------------------------
C  L O C A L  V A R I B L E S
C----------------------------------------------------------------
      INTEGER    I,J,K,IJK,IFAST(MVSIZ),NPRIN,INDEX(MVSIZ),NFAST
      INTEGER :: IFAST1
      INTEGER, DIMENSION(MVSIZ) :: IFAST1_INDEX
      my_real
     .    sv(3),ev(mvsiz,3),dirprv(mvsiz,3,3),ekk,dav,p3,g2(mvsiz),
     .    espin1(mvsiz), espin2(mvsiz), espin3(mvsiz), 
     .    espin4(mvsiz), espin5(mvsiz), espin6(mvsiz),
     .    e1, e2, e3, e4, e5, e6,espin(mvsiz,6), 
     .    e1q,e2q,e3q,e4q,e5q,e6q,
     .    e1t,e2t,e3t,e4t,e5t,e6t,b2(6),emin,emax
C----------------------------------------------------------------
         g2(1:nel)=two*g(1:nel)
          nprin = 0
          ifast1 = 0
          DO i=1,nel
            ifast(i) =1
            emin=min(es1(i),es2(i),es3(i),es4(i),es5(i),es6(i))
            emax=max(es1(i),es2(i),es3(i),es4(i),es5(i),es6(i))
            IF (emin<-em01.OR.emax>em01) THEN
             nprin = nprin + 1
             index(nprin) = i
             ifast(i) =0
             espin1(nprin)=es1(i)
             espin2(nprin)=es2(i)
             espin3(nprin)=es3(i)
             espin4(nprin)=es4(i)
             espin6(nprin)=es6(i)
             espin5(nprin)=es5(i)
            ELSE
              ifast1 = ifast1 + 1
              ifast1_index(ifast1) = i
            END IF
          ENDDO
C
#include "vectorize.inc"
          DO ijk=1,ifast1
           i = ifast1_index(ijk)
           e1=es1(i)
           e2=es2(i)
           e3=es3(i)
           e4=es4(i)
           e6=es6(i)
           e5=es5(i)
           e1q = e1*e1
           e2q = e2*e2
           e3q = e3*e3
           e4q = e4*e4
           e6q = e6*e6
           e5q = e5*e5
           b2(1)=e1q+e4q+e6q
           b2(2)=e4q+e2q+e5q
           b2(3)=e6q+e5q+e3q
           b2(4)=e1*e4+e2*e4+e5*e6
           b2(6)=e1*e6+e4*e5+e3*e6
           b2(5)=e4*e6+e2*e5+e3*e5
           es1(i)=half*(e1-half*b2(1))
           es2(i)=half*(e2-half*b2(2))
           es3(i)=half*(e3-half*b2(3))
           es4(i)=e4-half*b2(4)
           es6(i)=e6-half*b2(6)
           es5(i)=e5-half*b2(5)
          ENDDO
#include "vectorize.inc"
          DO ijk=1,ifast1
           i = ifast1_index(ijk)
            ekk=es1(i)+es2(i)+es3(i)-epsth(i)
            dav=-third*(es1(i)+es2(i)+es3(i))
            p3 = c1(i)*ekk
            sig(i,1)=p3+g2(i)*(es1(i)+dav)
            sig(i,2)=p3+g2(i)*(es2(i)+dav)
            sig(i,3)=p3+g2(i)*(es3(i)+dav)
            sig(i,4)=g(i)*es4(i)
            sig(i,5)=g(i)*es5(i)
            sig(i,6)=g(i)*es6(i)
          ENDDO
C
       IF (nprin>0) THEN
        CALL princ_u(
     1                  nprin    , espin1 , espin2 ,espin3 ,espin4 ,
     2                  espin5   ,espin6  ,ev   , dirprv  )
#include "vectorize.inc"
        DO k=1,nprin
          i = index(k)
          ekk=ev(k,1)+ev(k,2)+ev(k,3)
          dav=-third*ekk
          p3 = c1(i)*(ekk-epsth(i))
          sv(1)=p3 +g2(i)*(ev(k,1)+dav)
          sv(2)=p3 +g2(i)*(ev(k,2)+dav)
          sv(3)=p3 +g2(i)*(ev(k,3)+dav)
           sig(i,1) =  dirprv(k,1,1)*dirprv(k,1,1)*sv(1)
     .               + dirprv(k,1,2)*dirprv(k,1,2)*sv(2)
     .               + dirprv(k,1,3)*dirprv(k,1,3)*sv(3)     
           sig(i,2) =  dirprv(k,2,2)*dirprv(k,2,2)*sv(2)
     .               + dirprv(k,2,3)*dirprv(k,2,3)*sv(3)
     .               + dirprv(k,2,1)*dirprv(k,2,1)*sv(1)     
           sig(i,3) =  dirprv(k,3,3)*dirprv(k,3,3)*sv(3)
     .               + dirprv(k,3,1)*dirprv(k,3,1)*sv(1)
     .               + dirprv(k,3,2)*dirprv(k,3,2)*sv(2)     
           sig(i,4) =  dirprv(k,1,1)*dirprv(k,2,1)*sv(1)
     .               + dirprv(k,1,2)*dirprv(k,2,2)*sv(2)
     .               + dirprv(k,1,3)*dirprv(k,2,3)*sv(3)     
           sig(i,5) =  dirprv(k,2,2)*dirprv(k,3,2)*sv(2)
     .               + dirprv(k,2,3)*dirprv(k,3,3)*sv(3)
     .               + dirprv(k,2,1)*dirprv(k,3,1)*sv(1)     
           sig(i,6) =  dirprv(k,3,3)*dirprv(k,1,3)*sv(3)
     .               + dirprv(k,3,1)*dirprv(k,1,1)*sv(1)
     .               + dirprv(k,3,2)*dirprv(k,1,2)*sv(2)             
        ENDDO
       END IF
C      
      RETURN