844
845
846
847#include "implicit_f.inc"
848
849
850
851#include "mvsiz_p.inc"
852
853
854
855 INTEGER LLT
856
858 . dxdr(mvsiz), dydr(mvsiz), dzdr(mvsiz),
859 . dxdt(mvsiz), dydt(mvsiz), dzdt(mvsiz),
860 . drdx(mvsiz), dsdx(mvsiz), dtdx(mvsiz),
861 . drdy(mvsiz), dsdy(mvsiz), dtdy(mvsiz),
862 . drdz(mvsiz), dsdz(mvsiz), dtdz(mvsiz),
863 . xx(mvsiz,17) ,yy(mvsiz,17),zz(mvsiz,17),
864 . ni(mvsiz,17) ,rr(mvsiz) ,ss(mvsiz) ,tt(mvsiz)
865
866
867
868 INTEGER I,N
870 . dxds(mvsiz), dyds(mvsiz), dzds(mvsiz),
871 . dnidr(16),dnids(16),dnidt(16),
872 . d, aa, bb, det(mvsiz),r9 ,r13 ,s9 ,s10 ,s11 ,s12 ,t10 ,t14
874 . u_m_r,u_p_r,u_m_s,u_p_s,u_m_t,u_p_t,
875 . ums_umt,ums_upt,ups_umt,ups_upt,
876 . umr_ums,umr_ups,upr_ums,upr_ups,
877 . umt_umr,umt_upr,upt_umr,upt_upr,
878 . a,r05,s05,t05
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941 DO i=1,llt
942 r05 = half*rr(i)
943 s05 = half*ss(i)
944 t05 = half*tt(i)
945
946 u_m_r = half - r05
947 u_p_r = half + r05
948
949 u_m_s = half - s05
950 u_p_s = half + s05
951
952 u_m_t = half - t05
953 u_p_t = half + t05
954
955 ums_umt = u_m_s * u_m_t
956 ums_upt = u_m_s * u_p_t
957 ups_umt = u_p_s * u_m_t
958 ups_upt = u_p_s * u_p_t
959
960 umr_ums = u_m_r * u_m_s
961 umr_ups = u_m_r * u_p_s
962 upr_ums = u_p_r * u_m_s
963 upr_ups = u_p_r * u_p_s
964
965 umt_umr = u_m_t * u_m_r
966 umt_upr = u_m_t * u_p_r
967 upt_umr = u_p_t * u_m_r
968 upt_upr = u_p_t * u_p_r
969
970 a = -rr(i)-tt(i)-one
971 ni(i,1) = u_m_r * ums_umt * a
972 ni(i,5) = u_m_r * ups_umt * a
973 a = -rr(i)+tt(i)-one
974 ni(i,2) = u_m_r * ums_upt * a
975 ni(i,6) = u_m_r * ups_upt * a
976 a = rr(i)+tt(i)-one
977 ni(i,3) = u_p_r * ums_upt * a
978 ni(i,7) = u_p_r * ups_upt * a
979 a = rr(i)-tt(i)-one
980 ni(i,4) = u_p_r * ums_umt * a
981 ni(i,8) = u_p_r * ups_umt * a
982
983 a = -t05 - rr(i)
984 dnidr(1) = -ums_umt * a
985 dnidr(5) = -ups_umt * a
986 a = t05 - rr(i)
987 dnidr(2) = -ums_upt * a
988 dnidr(6) = -ups_upt * a
989 a = t05 + rr(i)
990 dnidr(3) = ums_upt * a
991 dnidr(7) = ups_upt * a
992 a = -t05 + rr(i)
993 dnidr(4) = ums_umt * a
994 dnidr(8) = ups_umt * a
995
996 a = -r05 - u_p_t
997 dnids(1) = -umt_umr * a
998 dnids(5) = umt_umr * a
999 a = -r05 - u_m_t
1000 dnids(2) = -upt_umr * a
1001 dnids(6) = upt_umr * a
1002 a = r05 - u_m_t
1003 dnids(3) = -upt_upr * a
1004 dnids(7) = upt_upr * a
1005 a = r05 - u_p_t
1006 dnids(4) = -umt_upr * a
1007 dnids(8) = umt_upr * a
1008
1009 a = -r05 - tt(i)
1010 dnidt(1) = -umr_ums * a
1011 dnidt(5) = -umr_ups * a
1012 a = -r05 + tt(i)
1013 dnidt(2) = umr_ums * a
1014 dnidt(6) = umr_ups * a
1015 a = +r05 + tt(i)
1016 dnidt(3) = upr_ums * a
1017 dnidt(7) = upr_ups * a
1018 a = +r05 - tt(i)
1019 dnidt(4) = -upr_ums * a
1020 dnidt(8) = -upr_ups * a
1021
1022 a = (one-rr(i)*rr(i))
1023 ni(i,10) = a * ums_upt
1024 ni(i,12) = a * ums_umt
1025 ni(i,14) = a * ups_upt
1026 ni(i,16) = a * ups_umt
1027
1028 a = half*a
1029 dnidt(10) = a * u_m_s
1030 dnidt(14) = a * u_p_s
1031 dnids(10) = -a * u_p_t
1032 dnids(12) = -a * u_m_t
1033
1034 a = -two*rr(i)
1035 dnidr(10) = a * ums_upt
1036 dnidr(12) = a * ums_umt
1037 dnidr(14) = a * ups_upt
1038 dnidr(16) = a * ups_umt
1039
1040 a = (one-tt(i)*tt(i))
1041 ni(i,9) = a * umr_ums
1042 ni(i,11) = a * upr_ums
1043 ni(i,13) = a * umr_ups
1044 ni(i,15) = a * upr_ups
1045
1046 ni(i,17) = -one
1047
1048 a = half*a
1049 dnidr(9) = -a * u_m_s
1050 dnidr(13) = -a * u_p_s
1051
1052 dnids(9) = -a * u_m_r
1053 dnids(11) = -a * u_p_r
1054
1055 a = -two*tt(i)
1056 dnidt(9) = a * umr_ums
1057 dnidt(11) = a * upr_ums
1058 dnidt(13) = a * umr_ups
1059 dnidt(15) = a * upr_ups
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069 dxdr(i) = dnidr(1)*xx(i,1) + dnidr(2)*xx(i,2) + dnidr(3)*xx(i,3)
1070 + + dnidr(4)*xx(i,4) + dnidr(5)*xx(i,5) + dnidr(6)*xx(i,6)
1071 + + dnidr(7)*xx(i,7) + dnidr(8)*xx(i,8)
1072 + + dnidr(9)*(xx(i,9) - xx(i,11)) + dnidr(10)*xx(i,10)
1073 + + dnidr(12)*xx(i,12) + dnidr(13)*(xx(i,13) - xx(i,15))
1074 + + dnidr(14)*xx(i,14) + dnidr(16)*xx(i,16)
1075
1076 dxds(i) = dnids(1)*xx(i,1) + dnids(2)*xx(i,2) + dnids(3)*xx(i,3)
1077 + + dnids(4)*xx(i,4) + dnids(5)*xx(i,5) + dnids(6)*xx(i,6)
1078 + + dnids(7)*xx(i,7) + dnids(8)*xx(i,8)
1079 + + dnids(9)* (xx(i,9) - xx(i,13))
1080 + + dnids(10)*(xx(i,10) - xx(i,14))
1081 + + dnids(11)*(xx(i,11) - xx(i,15))
1082 + + dnids(12)*(xx(i,12) - xx(i,16))
1083
1084 dxdt(i) = dnidt(1)*xx(i,1) + dnidt(2)*xx(i,2) + dnidt(3)*xx(i,3)
1085 + + dnidt(4)*xx(i,4) + dnidt(5)*xx(i,5) + dnidt(6)*xx(i,6)
1086 + + dnidt(7)*xx(i,7) + dnidt(8)*xx(i,8)
1087 + + dnidt(9)*xx(i,9) + dnidt(10)*(xx(i,10) - xx(i,12))
1088 + + dnidt(11)*xx(i,11) + dnidt(13)*xx(i,13)
1089 + + dnidt(14)*(xx(i,14) - xx(i,16)) + dnidt(15)*xx(i,15)
1090
1091
1092
1093 dydr(i) = dnidr(1)*yy(i,1) + dnidr(2)*yy(i,2) + dnidr(3)*yy(i,3)
1094 + + dnidr(4)*yy(i,4) + dnidr(5)*yy(i,5) + dnidr(6)*yy(i,6)
1095 + + dnidr(7)*yy(i,7) + dnidr(8)*yy(i,8)
1096 + + dnidr(9)*(yy(i,9) - yy(i,11)) + dnidr(10)*yy(i,10)
1097 + + dnidr(12)*yy(i,12) + dnidr(13)*(yy(i,13) - yy(i,15))
1098 + + dnidr(14)*yy(i,14) + dnidr(16)*yy(i,16)
1099
1100 dyds(i) = dnids(1)*yy(i,1) + dnids(2)*yy(i,2) + dnids(3)*yy(i,3)
1101 + + dnids(4)*yy(i,4) + dnids(5)*yy(i,5) + dnids(6)*yy(i,6)
1102 + + dnids(7)*yy(i,7) + dnids(8)*yy(i,8)
1103 + + dnids(9)* (yy(i,9) - yy(i,13))
1104 + + dnids(10)*(yy(i,10) - yy(i,14))
1105 + + dnids(11)*(yy(i,11) - yy(i,15))
1106 + + dnids(12)*(yy(i,12) - yy(i,16))
1107
1108 dydt(i) = dnidt(1)*yy(i,1) + dnidt(2)*yy(i,2) + dnidt(3)*yy(i,3)
1109 + + dnidt(4)*yy(i,4) + dnidt(5)*yy(i,5) + dnidt(6)*yy(i,6)
1110 + + dnidt(7)*yy(i,7) + dnidt(8)*yy(i,8)
1111 + + dnidt(9)*yy(i,9) + dnidt(10)*(yy(i,10) - yy(i,12))
1112 + + dnidt(11)*yy(i,11) + dnidt(13)*yy(i,13)
1113 + + dnidt(14)*(yy(i,14) - yy(i,16)) + dnidt(15)*yy(i,15)
1114
1115
1116
1117 dzdr(i) = dnidr(1)*zz(i,1) + dnidr(2)*zz(i,2) + dnidr(3)*zz(i,3)
1118 + + dnidr(4)*zz(i,4) + dnidr(5)*zz(i,5) + dnidr(6)*zz(i,6)
1119 + + dnidr(7)*zz(i,7) + dnidr(8)*zz(i,8)
1120 + + dnidr(9)*(zz(i,9) - zz(i,11)) + dnidr(10)*zz(i,10)
1121 + + dnidr(12)*zz(i,12) + dnidr(13)*(zz(i,13) - zz(i,15))
1122 + + dnidr(14)*zz(i,14) + dnidr(16)*zz(i,16)
1123
1124 dzds(i) = dnids(1)*zz(i,1) + dnids(2)*zz(i,2) + dnids(3)*zz(i,3)
1125 + + dnids(4)*zz(i,4) + dnids(5)*zz(i,5) + dnids(6)*zz(i,6)
1126 + + dnids(7)*zz(i,7) + dnids(8)*zz(i,8)
1127 + + dnids(9)* (zz(i,9) - zz(i,13))
1128 + + dnids(10)*(zz(i,10) - zz(i,14))
1129 + + dnids(11)*(zz(i,11) - zz(i,15))
1130 + + dnids(12)*(zz(i,12) - zz(i,16))
1131
1132 dzdt(i) = dnidt(1)*zz(i,1) + dnidt(2)*zz(i,2) + dnidt(3)*zz(i,3)
1133 + + dnidt(4)*zz(i,4) + dnidt(5)*zz(i,5) + dnidt(6)*zz(i,6)
1134 + + dnidt(7)*zz(i,7) + dnidt(8)*zz(i,8)
1135 + + dnidt(9)*zz(i,9) + dnidt(10)*(zz(i,10) - zz(i,12))
1136 + + dnidt(11)*zz(i,11) + dnidt(13)*zz(i,13)
1137 + + dnidt(14)*(zz(i,14) - zz(i,16)) + dnidt(15)*zz(i,15)
1138
1139
1140
1141
1142 drdx(i)=dyds(i)*dzdt(i)-dzds(i)*dydt(i)
1143 drdy(i)=dzds(i)*dxdt(i)-dxds(i)*dzdt(i)
1144 drdz(i)=dxds(i)*dydt(i)-dyds(i)*dxdt(i)
1145
1146 dsdz(i)=dxdt(i)*dydr(i)-dydt(i)*dxdr(i)
1147 dsdy(i)=dzdt(i)*dxdr(i)-dxdt(i)*dzdr(i)
1148 dsdx(i)=dydt(i)*dzdr(i)-dzdt(i)*dydr(i)
1149
1150 dtdx(i)=dydr(i)*dzds(i)-dzdr(i)*dyds(i)
1151 dtdy(i)=dzdr(i)*dxds(i)-dxdr(i)*dzds(i)
1152 dtdz(i)=dxdr(i)*dyds(i)-dydr(i)*dxds(i)
1153
1154 det(i) = dxdr(i) * drdx(i)
1155 . + dydr(i) * drdy(i)
1156 . + dzdr(i) * drdz(i)
1157
1158
1159
1160
1161
1162
1163
1164 ENDDO
1165
1166 DO i=1,llt
1167
1168
1169
1170
1171 d = one/det(i)
1172 drdx(i)=d*drdx(i)
1173 dsdx(i)=d*dsdx(i)
1174 dtdx(i)=d*dtdx(i)
1175
1176 drdy(i)=d*drdy(i)
1177 dsdy(i)=d*dsdy(i)
1178 dtdy(i)=d*dtdy(i)
1179
1180 drdz(i)=d*drdz(i)
1181 dsdz(i)=d*dsdz(i)
1182 dtdz(i)=d*dtdz(i)
1183
1184
1185
1186
1187
1188
1189 ENDDO
1190
1191 RETURN