Esope source files

   1 : C FFFACE    SOURCE    CHAT      05/01/12    23:58:45     5004           
   2 :       SUBROUTINE FFFACE(XA,NPF)
   3 :       IMPLICIT INTEGER(I-N)
   4 :       IMPLICIT REAL*8 (A-H,O-Z)                                                 
   5 : C************************************************************************       
   6 : C                                                                               
   7 : C                                                                               
   8 : C                                                                               
   9 : C************************************************************************       
  10 :       DIMENSION FN(20),XA(3,21),COOR(3,7)
  11 : C                                                                               
  12 :       IF(NPF.EQ.8)THEN
  13 :       X=0.5D0
  14 :       Y=0.5D0
  15 : 
  16 :       FN(1)=-2.D0*(1.D0-X)*(1.D0-Y)*(X+Y-0.5D0)
  17 :       FN(2)= 4.D0*X*(1.D0-X)*(1.D0-Y)
  18 :       FN(3)=-2.D0*X*(1.D0-Y)*(Y-X+0.5D0)
  19 :       FN(4)= 4.D0*X*Y*(1.D0-Y)
  20 :       FN(5)=-2.D0*X*Y*(1.5D0-X-Y)
  21 :       FN(6)= 4.D0*X*Y*(1.D0-X)
  22 :       FN(7)=-2.D0*Y*(1.D0-X)*(X-Y+0.5D0)
  23 :       FN(8)= 4.D0*Y*(1.D0-Y)*(1.D0-X)
  24 : 
  25 :       DO 2 N=1,3
  26 :       XA(N,NPF+1)=
  27 :      &         FN(1)*XA(N,1)+FN(2)*XA(N,2)+FN(3)*XA(N,3)+FN(4)*XA(N,4)
  28 :      &+        FN(5)*XA(N,5)+FN(6)*XA(N,6)+FN(7)*XA(N,7)+FN(8)*XA(N,8)
  29 :  2    CONTINUE
  30 : 
  31 :       ELSEIF(NPF.EQ.6)THEN
  32 :       A=SQRT(2.D0)
  33 :       AS2=A/2.D0
  34 :       X=AS2
  35 :       Y=AS2
  36 :       FN(1)=(X+Y-AS2)*(X+Y-A)
  37 :       FN(2)=-2.D0*X*(X+Y-A)
  38 :       FN(3)=-X*(X-AS2)
  39 :       FN(4)=2.D0*X*Y
  40 :       FN(5)=Y*(Y-AS2)
  41 :       FN(6)=-2.D0*Y*(X+Y-A)
  42 : 
  43 :       DO 4 N=1,3
  44 :       XA(N,NPF+1)=(XA(N,1)+XA(N,3)+XA(N,5))/3.D0
  45 : C?    XA(N,NPF+1)=
  46 : C?   &         FN(1)*XA(N,1)+FN(2)*XA(N,2)+FN(3)*XA(N,3)+FN(4)*XA(N,4)
  47 : C?   &+        FN(5)*XA(N,5)+FN(6)*XA(N,6)
  48 :  4    CONTINUE
  49 : 
  50 :       ELSEIF(NPF.EQ.20)THEN
  51 : 
  52 :       X=0.5D0
  53 :       Y=0.5D0
  54 :       Z=0.5D0
  55 : 
  56 :       FN(1)=-2.D0*(1.D0-X)*(1.D0-Y)*(1.D0-Z)*(X+Y+Z-0.5D0)
  57 :       FN(2)= 4.D0*X*(1.D0-X)*(1.D0-Y)*(1.D0-Z)
  58 :       FN(3)=-2.D0*X*(1.D0-Y)*(1.D0-Z)*(Z+Y-X+0.5D0)
  59 :       FN(4)= 4.D0*X*Y*(1.D0-Y)*(1.D0-Z)
  60 :       FN(5)=-2.D0*X*Y*(1.D0-Z)*(1.5D0-X-Y+Z)
  61 :       FN(6)= 4.D0*X*Y*(1.D0-X)*(1.D0-Z)
  62 :       FN(7)=-2.D0*Y*(1.D0-X)*(1.D0-Z)*(X-Y+Z+0.5D0)
  63 :       FN(8)= 4.D0*Y*(1.D0-Y)*(1.D0-X)*(1.D0-Z)
  64 :       FN(9) = 4.D0*Z*(1.D0-Z)*(1.D0-X)*(1.D0-Y)
  65 :       FN(10)= 4.D0*X*Z*(1.D0-Z)*(1.D0-Y)
  66 :       FN(11)= 4.D0*X*Y*Z*(1.D0-Z)
  67 :       FN(12)= 4.D0*Y*Z*(1.D0-Z)*(1.D0-X)
  68 :       FN(13)=-2.D0*(1.D0-X)*(1.D0-Y)*Z*(X+Y-Z+0.5D0)
  69 :       FN(14)= 4.D0*X*(1.D0-X)*(1.D0-Y)*Z
  70 :       FN(15)=-2.D0*X*(1.D0-Y)*Z*(-Z+Y-X+1.5D0)
  71 :       FN(16)= 4.D0*X*Y*(1.D0-Y)*Z
  72 :       FN(17)=-2.D0*X*Y*Z*(2.5D0-X-Y-Z)
  73 :       FN(18)= 4.D0*X*Y*(1.D0-X)*Z
  74 :       FN(19)=-2.D0*Y*(1.D0-X)*Z*(X-Y-Z+1.5D0)
  75 :       FN(20)= 4.D0*Y*(1.D0-Y)*(1.D0-X)*Z
  76 : 
  77 :       DO 20 N=1,3
  78 :       XA(N,NPF+1)=
  79 :      &         FN(1)*XA(N,1)+FN(2)*XA(N,2)+FN(3)*XA(N,3)+FN(4)*XA(N,4)
  80 :      &+        FN(5)*XA(N,5)+FN(6)*XA(N,6)+FN(7)*XA(N,7)+FN(8)*XA(N,8)
  81 :      &+  FN(9)*XA(N,9)+FN(10)*XA(N,10)+FN(11)*XA(N,11)+FN(12)*XA(N,12)
  82 :      &+FN(13)*XA(N,13)+FN(14)*XA(N,14)+FN(15)*XA(N,15)+FN(16)*XA(N,16)
  83 :      &+FN(17)*XA(N,17)+FN(18)*XA(N,18)+FN(19)*XA(N,19)+FN(20)*XA(N,20)
  84 :  20   CONTINUE
  85 : 
  86 :       ELSEIF(NPF.EQ.15)THEN
  87 : 
  88 :       A=SQRT(2.D0)
  89 :       AS2=A*0.5D0
  90 :       CA2=1.D0
  91 :       CA4=2.D0
  92 :       A4=4.D0/A
  93 : 
  94 : C     COOR(1,1)=AS2
  95 : C     COOR(2,1)=0.D0
  96 : C     COOR(3,1)=0.5D0
  97 : 
  98 :       COOR(1,2)=AS2
  99 :       COOR(2,2)=AS2
 100 :       COOR(3,2)=0.5D0
 101 : 
 102 : C     COOR(1,3)=0.D0
 103 : C     COOR(2,3)=AS2
 104 : C     COOR(3,3)=0.5D0
 105 : 
 106 :       X=COOR(1,2)
 107 :       Y=COOR(2,2)
 108 :       Z=COOR(3,2)
 109 : C     write(6,*)' X,Y,Z=',X,Y,Z
 110 : 
 111 :       FN(1)= CA2*(X+Y-A)*(1.D0-Z)*(X+Y+A*Z-AS2)
 112 :       FN(2)=-CA4*X*(X+Y-A)*(1.D0-Z)
 113 :       FN(3)=-CA2*X*(1.D0-Z)*(-X+A*Z+AS2)
 114 :       FN(4)= CA4*X*Y*(1.D0-Z)
 115 :       FN(5)=-CA2*Y*(1.D0-Z)*(-Y+A*Z+AS2)
 116 :       FN(6)=-CA4*Y*(X+Y-A)*(1.D0-Z)
 117 :       FN(7)=-A4*Z*(1.D0-Z)*(X+Y-A)
 118 :       FN(8)= A4*X*Z*(1.D0-Z)
 119 :       FN(9) = A4*Y*Z*(1.D0-Z)
 120 :       FN(10)= CA2*(X+Y-A)*Z*(X+Y+A*(1.D0-Z)-AS2)
 121 :       FN(11)=-CA4*X*(X+Y-A)*Z
 122 :       FN(12)=-CA2*X*Z*(-X+A*(1.D0-Z)+AS2)
 123 :       FN(13)= CA4*X*Y*Z
 124 :       FN(14)=-CA2*Y*Z*(-Y+A*(1.D0-Z)+AS2)
 125 :       FN(15)=-CA4*Y*Z*(X+Y-A)
 126 : 
 127 :       DO 15 N=1,3
 128 :       XA(N,NPF+1)=(XA(N,7)+XA(N,8)+XA(N,9))/3.D0
 129 : C?    XA(N,NPF+1)=
 130 : C?   &         FN(1)*XA(N,1)+FN(2)*XA(N,2)+FN(3)*XA(N,3)+FN(4)*XA(N,4)
 131 : C?   &+        FN(5)*XA(N,5)+FN(6)*XA(N,6)+FN(7)*XA(N,7)+FN(8)*XA(N,8)
 132 : C?   &+  FN(9)*XA(N,9)+FN(10)*XA(N,10)+FN(11)*XA(N,11)+FN(12)*XA(N,12)
 133 : C?   &+FN(13)*XA(N,13)+FN(14)*XA(N,14)+FN(15)*XA(N,15)
 134 :  15   CONTINUE
 135 :       ELSE
 136 : 
 137 :       DO 30 N=1,3
 138 :       U=0.D0
 139 :       DO 31 I=1,NPF
 140 :       U=U+XA(N,I)
 141 :  31   CONTINUE
 142 :       XA(N,NPF+1)=U/FLOAT(NPF)
 143 :  30   CONTINUE
 144 : 
 145 :       ENDIF
 146 : 
 147 : 
 148 : 
 149 : 1002  format(10(1x,1pe11.4))
 150 :       RETURN                                                                    
 151 :       END                                                                       
 152 :  
 153 :