1043
1044
1045
1046#include "implicit_f.inc"
1047
1048
1049
1050#include "mvsiz_p.inc"
1051
1052
1053
1054 INTEGER LLT
1055
1057 . dxdr(mvsiz), dydr(mvsiz), dzdr(mvsiz),
1058 . dxds(mvsiz), dyds(mvsiz), dzds(mvsiz),
1059 . dxdt(mvsiz), dydt(mvsiz), dzdt(mvsiz),
1060 . drdx(mvsiz), dsdx(mvsiz), dtdx(mvsiz),
1061 . drdy(mvsiz), dsdy(mvsiz), dtdy(mvsiz),
1062 . drdz(mvsiz), dsdz(mvsiz), dtdz(mvsiz),
1063 . xx(mvsiz,21) ,yy(mvsiz,21),zz(mvsiz,21),
1064 . ni(mvsiz,21) ,rr(mvsiz) ,ss(mvsiz) ,tt(mvsiz)
1065
1066
1067
1068 INTEGER I,N
1070 . dnidr(20),dnids(20),dnidt(20),
1071 . d, aa, bb, det(mvsiz),r9 ,r13 ,s9 ,s10 ,s11 ,s12 ,t10 ,t14
1073 . u_m_r,u_p_r,u_m_s,u_p_s,u_m_t,u_p_t,
1074 . ums_umt,ums_upt,ups_umt,ups_upt,
1075 . umr_ums,umr_ups,upr_ums,upr_ups,
1076 . umt_umr,umt_upr,upt_umr,upt_upr,
1077 . a,r05,s05,t05
1078
1079
1080
1081
1082
1083
1084
1085
1086
1087
1088
1089
1090
1091
1092
1093
1094
1095
1096
1097
1098
1099
1100
1101
1102
1103
1104
1105
1106
1107
1108
1109
1110
1111
1112
1113
1114
1115
1116
1117
1118
1119
1120
1121
1122
1123
1124
1125
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
1149 DO i=1,llt
1150 r05 = half*rr(i)
1151 s05 = half*ss(i)
1152 t05 = half*tt(i)
1153
1154 u_m_r = half - r05
1155 u_p_r = half + r05
1156
1157 u_m_s = half - s05
1158 u_p_s = half + s05
1159
1160 u_m_t = half - t05
1161 u_p_t = half + t05
1162
1163 ums_umt = u_m_s * u_m_t
1164 ums_upt = u_m_s * u_p_t
1165 ups_umt = u_p_s * u_m_t
1166 ups_upt = u_p_s * u_p_t
1167
1168 umr_ums = u_m_r * u_m_s
1169 umr_ups = u_m_r * u_p_s
1170 upr_ums = u_p_r * u_m_s
1171 upr_ups = u_p_r * u_p_s
1172
1173 umt_umr = u_m_t * u_m_r
1174 umt_upr = u_m_t * u_p_r
1175 upt_umr = u_p_t * u_m_r
1176 upt_upr = u_p_t * u_p_r
1177
1178
1179
1180
1181 ni(i,1) = u_m_r * ums_umt * (-rr(i)-ss(i)-tt(i
1182 ni(i,2) = u_m_r * ums_upt * (-rr(i)-ss(i)+tt(i)-two)
1183 ni(i,3) = u_p_r * ums_upt * ( rr(i)-ss(i)+tt(i)-two)
1184 ni(i,4) = u_p_r * ums_umt * ( rr(i)-ss(i)-tt(i)-two)
1185 ni(i,5) = u_m_r * ups_umt * (-rr(i)+ss(i)-tt(i)-two)
1186 ni(i,6) = u_m_r * ups_upt * (-rr(i)+ss(i)+tt(i)-two)
1187 ni(i,7) = u_p_r * ups_upt * ( rr(i)+ss(i)+tt(i)-two)
1188 ni(i,8) = u_p_r * ups_umt * ( rr(i)+ss(i)-tt(i)-two)
1189
1190 dnidr(1) = -ums_umt * (u_m_s + u_m_t - rr(i) -three_half)
1191 dnidr(2) = -ums_upt * (u_m_s + u_p_t - rr(i) -three_half)
1192 dnidr(3) = ums_upt * (u_m_s + u_p_t + rr(i) -three_half)
1193 dnidr(4) = ums_umt * (u_m_s + u_m_t + rr(i) -three_half)
1194 dnidr(5) = -ups_umt * (u_p_s + u_m_t - rr(i) -three_half)
1195 dnidr(6) = -ups_upt * (u_p_s + u_p_t - rr(i) -three_half)
1196 dnidr(7) = ups_upt * (u_p_s + u_p_t + rr(i) -three_half)
1197 dnidr(8) = ups_umt * (u_p_s + u_m_t + rr(i) -three_half)
1198
1199 dnids(1) = -umt_umr * (u_m_r + u_m_t - ss(i) -three_half)
1200 dnids(2) = -upt_umr * (u_m_r + u_p_t - ss(i) -three_half)
1201 dnids(3) = -upt_upr * (u_p_r + u_p_t - ss(i) -three_half)
1202 dnids(4) = -umt_upr * (u_p_r + u_m_t - ss(i) -three_half)
1203 dnids(5) = umt_umr * (u_m_r + u_m_t + ss(i) -three_half)
1204 dnids(6) = upt_umr * (u_m_r + u_p_t + ss(i) -three_half)
1205 dnids(7) = upt_upr * (u_p_r + u_p_t + ss(i) -three_half)
1206 dnids(8) = umt_upr * (u_p_r + u_m_t + ss(i) -three_half)
1207
1208 dnidt(1) = -umr_ums * (u_m_r + u_m_s - tt(i) -three_half)
1209 dnidt(2) = umr_ums * (u_m_r + u_m_s + tt(i) -three_half)
1210 dnidt(3) = upr_ums * (u_p_r + u_m_s + tt(i) -three_half)
1211 dnidt(4) = -upr_ums * (u_p_r + u_m_s - tt(i) -three_half)
1212 dnidt(5) = -umr_ups * (u_m_r + u_p_s - tt(i) -three_half)
1213 dnidt(6) = umr_ups * (u_m_r + u_p_s + tt(i) -three_half)
1214 dnidt(7) = upr_ups * (u_p_r + u_p_s + tt(i) -three_half)
1215 dnidt(8) = -upr_ups * (u_p_r + u_p_s - tt(i) -three_half)
1216
1217 a = (one-rr(i)*rr(i))
1218 ni(i,10) = a * ums_upt
1219 ni(i,12) = a * ums_umt
1220 ni(i,18) = a * ups_upt
1221 ni(i,20) = a * ups_umt
1222
1223 a = half*a
1224 dnidt(10) = a * u_m_s
1225 dnidt(18) = a * u_p_s
1226 dnidt(12) = -dnidt(10)
1227 dnidt(20) = -dnidt(18)
1228
1229 dnids(18) = a * u_p_t
1230 dnids(20) = a * u_m_t
1231 dnids(10) = -dnids(18)
1232 dnids(12) = -dnids(20)
1233
1234 a = -two*rr(i)
1235 dnidr(10) = a * ums_upt
1236 dnidr(12) = a * ums_umt
1237 dnidr(18) = a * ups_upt
1238 dnidr(20) = a * ups_umt
1239
1240 a = (one-ss(i)*ss(i))
1241 ni(i,13) = a * umt_umr
1242 ni(i,14) = a * upt_umr
1243 ni(i,15) = a * upt_upr
1244 ni(i,16) = a * umt_upr
1245
1246 a = half*a
1247 dnidr(15) = a * u_p_t
1248 dnidr(16) = a * u_m_t
1249 dnidr(13) = -dnidr(16)
1250 dnidr(14) = -dnidr(15)
1251
1252 dnidt(14) = a * u_m_r
1253 dnidt(15) = a * u_p_r
1254 dnidt(13) = -dnidt(14)
1255 dnidt(16) = -dnidt(15)
1256
1257 a = -two*ss(i)
1258 dnids(13) = a * umt_umr
1259 dnids(14) = a * upt_umr
1260 dnids(15) = a * upt_upr
1261 dnids(16) = a * umt_upr
1262
1263 a = (one-tt(i)*tt(i))
1264 ni(i,9) = a * umr_ums
1265 ni(i,11) = a * upr_ums
1266 ni(i,17) = a * umr_ups
1267 ni(i,19) = a * upr_ups
1268
1269 ni(i,21) = -one
1270
1271 a = half*a
1272 dnidr(11) = a * u_m_s
1273 dnidr(19) = a * u_p_s
1274 dnidr(9) = -dnidr(11)
1275 dnidr(17) = -dnidr(19)
1276
1277 dnids(17) = a * u_m_r
1278 dnids(19) = a * u_p_r
1279 dnids(9) = -dnids(17)
1280 dnids(11) = -dnids(19)
1281
1282 a = -two*tt(i)
1283 dnidt(9) = a * umr_ums
1284 dnidt(11) = a * upr_ums
1285 dnidt(17) = a * umr_ups
1286 dnidt(19) = a * upr_ups
1287
1288
1289
1290
1291
1292
1293
1294
1295
1296
1297 dxdr(i) = dnidr(1)*xx(i,1) + dnidr(2)*xx(i,2) + dnidr(3)*xx(i,3)
1298 + + dnidr(4)*xx(i,4) + dnidr(5)*xx(i,5) + dnidr(6)*xx(i,6)
1299 + + dnidr(7)*xx(i,7) + dnidr(8)*xx(i,8)
1300 + + dnidr(9)*xx(i,9) + dnidr(10)*xx(i,10) + dnidr(11)*xx(i,11)
1301 + + dnidr(12)*xx(i,12) + dnidr(13)*xx(i,13) + dnidr(14)*xx(i,14)
1302 + + dnidr(15)*xx(i,15) + dnidr(16)*xx(i,16) + dnidr(17)*xx(i,17)
1303 + + dnidr(18)*xx(i,18) + dnidr(19)*xx(i,19) + dnidr(20)*xx(i,20)
1304
1305 dxds(i) = dnids(1)*xx(i,1) + dnids(2)*xx(i,2) + dnids(3)*xx(i,3)
1306 + + dnids(4)*xx(i,4) + dnids(5)*xx(i,5) + dnids(6)*xx(i,6)
1307 + + dnids(7)*xx(i,7) + dnids(8)*xx(i,8)
1308 + + dnids(9)*xx(i,9) + dnids(10)*xx(i,10) + dnids(11)*xx(i,11)
1309 + + dnids(12)*xx(i,12) + dnids(13)*xx(i,13) + dnids(14)*xx(i,14)
1310 + + dnids(15)*xx(i,15) + dnids(16)*xx(i,16) + dnids(17)*xx(i,17)
1311 + + dnids(18)*xx(i,18) + dnids(19)*xx(i,19) + dnids(20)*xx(i,20)
1312
1313 dxdt(i) = dnidt(1)*xx(i,1) + dnidt(2)*xx(i,2) + dnidt(3)*xx(i,3)
1314 + + dnidt(4)*xx(i,4) + dnidt(5)*xx(i,5) + dnidt(6)*xx(i,6)
1315 + + dnidt(7)*xx(i,7) + dnidt(8)*xx(i,8)
1316 + + dnidt(9)*xx(i,9) + dnidt(10)*xx(i,10) + dnidt(11)*xx(i,11)
1317 + + dnidt(12)*xx(i,12) + dnidt(13)*xx(i,13) + dnidt(14)*xx(i,14)
1318 + + dnidt(15)*xx(i,15) + dnidt(16)*xx(i,16) + dnidt(17)*xx(i,17)
1319 + + dnidt(18)*xx(i,18) + dnidt(19)*xx(i,19) + dnidt(20)*xx(i,20)
1320
1321
1322
1323 dydr(i) = dnidr(1)*yy(i,1) + dnidr(2)*yy(i,2) + dnidr(3)*yy(i,3)
1324 + + dnidr(4)*yy(i,4) + dnidr(5)*yy(i,5) + dnidr(6)*yy(i,6)
1325 + + dnidr(7)*yy(i,7) + dnidr(8)*yy(i,8)
1326 + + dnidr(9)*yy(i,9) + dnidr(10)*yy(i,10) + dnidr(11)*yy(i,11)
1327 + + dnidr(12)*yy(i,12) + dnidr(13)*yy(i,13) + dnidr(14)*yy(i,14)
1328 + + dnidr(15)*yy(i,15) + dnidr(16)*yy(i,16) + dnidr(17)*yy(i,17)
1329 + + dnidr(18)*yy(i,18) + dnidr(19)*yy(i,19) + dnidr(20)*yy(i,20)
1330
1331 dyds(i) = dnids(1)*yy(i,1) + dnids(2)*yy(i,2) + dnids(3)*yy(i,3)
1332 + + dnids(4)*yy(i,4) + dnids(5)*yy(i,5) + dnids(6)*yy(i,6)
1333 + + dnids(7)*yy(i,7) + dnids(8)*yy(i,8)
1334 + + dnids(9)*yy(i,9) + dnids(10)*yy(i,10) + dnids(11)*yy(i,11)
1335 + + dnids(12)*yy(i,12) + dnids(13)*yy(i,13) + dnids(14)*yy(i,14)
1336 + + dnids(15)*yy(i,15) + dnids(16)*yy(i,16) + dnids(17)*yy(i,17)
1337 + + dnids(18)*yy(i,18) + dnids(19)*yy(i,19) + dnids(20)*yy(i,20)
1338
1339 dydt(i) = dnidt(1)*yy(i,1) + dnidt(2)*yy(i,2) + dnidt(3)*yy(i,3)
1340 + + dnidt(4)*yy(i,4) + dnidt(5)*yy(i,5) + dnidt(6)*yy(i,6)
1341 + + dnidt(7)*yy(i,7) + dnidt(8)*yy(i,8)
1342 + + dnidt(9)*yy(i,9) + dnidt(10)*yy(i,10) + dnidt(11)*yy(i,11)
1343 + + dnidt(12)*yy(i,12) + dnidt(13)*yy(i,13) + dnidt(14)*yy(i,14)
1344 + + dnidt(15)*yy(i,15) + dnidt(16)*yy(i,16) + dnidt(17)*yy(i,17)
1345 + + dnidt(18)*yy(i,18) + dnidt(19)*yy(i,19) + dnidt(20)*yy(i,20)
1346
1347
1348
1349 dzdr(i) = dnidr(1)*zz(i,1) + dnidr(2)*zz(i,2) + dnidr(3)*zz(i,3)
1350 + + dnidr(4)*zz(i,4) + dnidr(5)*zz(i,5) + dnidr(6)*zz(i,6)
1351 + + dnidr(7)*zz(i,7) + dnidr(8)*zz(i,8)
1352 + + dnidr(9)*zz(i,9) + dnidr(10)*zz(i,10) + dnidr(11)*zz(i,11)
1353 + + dnidr(12)*zz(i,12) + dnidr(13)*zz(i,13) + dnidr(14)*zz(i,14)
1354 + + dnidr(15)*zz(i,15) + dnidr(16)*zz(i,16) + dnidr(17)*zz(i,17)
1355 + + dnidr(18)*zz(i,18) + dnidr(19)*zz(i,19) + dnidr(20)*zz(i,20)
1356
1357 dzds(i) = dnids(1)*zz(i,1) + dnids(2)*zz(i,2) + dnids(3)*zz(i,3)
1358 + + dnids(4)*zz(i,4) + dnids(5)*zz(i,5) + dnids(6)*zz(i,6)
1359 + + dnids(7)*zz(i,7) + dnids(8)*zz(i,8)
1360 + + dnids(9)*zz(i,9) + dnids(10)*zz(i,10) + dnids(11)*zz(i,11)
1361 + + dnids(12)*zz(i,12) + dnids(13)*zz(i,13) + dnids(14)*zz(i,14)
1362 + + dnids(15)*zz(i,15) + dnids(16)*zz(i,16) + dnids(17)*zz(i,17)
1363 + + dnids(18)*zz(i,18) + dnids(19)*zz(i,19) + dnids(20)*zz(i,20)
1364
1365 dzdt(i) = dnidt(1)*zz(i,1) + dnidt(2)*zz(i,2) + dnidt(3)*zz(i,3)
1366 + + dnidt(4)*zz(i,4) + dnidt(5)*zz(i,5) + dnidt(6)*zz(i,6)
1367 + + dnidt(7)*zz(i,7) + dnidt(8)*zz(i,8)
1368 + + dnidt(9)*zz(i,9) + dnidt(10)*zz(i,10) + dnidt(11)*zz(i,11)
1369 + + dnidt(12)*zz(i,12) + dnidt(13)*zz(i,13) + dnidt(14)*zz(i,14)
1370 + + dnidt(15)*zz(i,15) + dnidt(16)*zz(i,16) + dnidt(17)*zz(i,17)
1371 + + dnidt(18)*zz(i,18) + dnidt(19)*zz(i,19) + dnidt(20)*zz(i,20)
1372
1373
1374
1375
1376 drdx(i)=dyds(i)*dzdt(i)-dzds(i)*dydt(i)
1377 drdy(i)=dzds(i)*dxdt(i)-dxds(i)*dzdt(i)
1378 drdz(i)=dxds(i)*dydt(i)-dyds(i)*dxdt(i)
1379
1380 dsdz(i)=dxdt(i)*dydr(i)-dydt(i)*dxdr(i)
1381 dsdy(i)=dzdt(i)*dxdr(i)-dxdt(i)*dzdr(i)
1382 dsdx(i)=dydt(i)*dzdr(i)-dzdt(i)*dydr(i)
1383
1384 dtdx(i)=dydr(i)*dzds(i)-dzdr(i)*dyds(i)
1385 dtdy(i)=dzdr(i)*dxds(i)-dxdr(i)*dzds(i)
1386 dtdz(i)=dxdr(i)*dyds(i)-dydr(i)*dxds(i)
1387
1388 det(i) = dxdr(i) * drdx(i)
1389 . + dydr(i) * drdy(i)
1390 . + dzdr(i) * drdz(i)
1391
1392
1393
1394
1395
1396
1397
1398 ENDDO
1399
1400 DO i=1,llt
1401
1402
1403
1404
1405 d = one/det(i)
1406 drdx(i)=d*drdx(i)
1407 dsdx(i)=d*dsdx(i)
1408 dtdx(i)=d*dtdx(i)
1409
1410 drdy(i)=d*drdy(i)
1411 dsdy(i)=d*dsdy(i)
1412 dtdy(i)=d*dtdy(i)
1413
1414 drdz(i)=d*drdz(i)
1415 dsdz(i)=d*dsdz(i)
1416 dtdz(i)=d*dtdz(i)
1417
1418
1419
1420
1421
1422
1423 ENDDO
1424
1425 RETURN