# 図7 s2 曲面と c2 曲線 TMP1="tmp1$$" TMP2="tmp2$$" trap 'rm -f $TMP1 $TMP2; exit 0' 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 /usr/bin/awk ' BEGIN { PI = atan2(0, -1) PI2 = atan2(0, -1) * 2 a = 25 b = 10 k1 = sqrt(1 / 1.0 + 1) k2 = sqrt(1 / 6^2 + 1) p0 = 0 N = 0 Ns = 5 Np = 20 pmin = p0 pmax = p0 + PI * 3 # s2 線 for (t = 0; t <= PI2; t += PI / Ns) { p1 = atan2(k1*sqrt(k1^2 - 1)*cos(t), sin(t)) p2 = p1 + PI while (p1 < pmax) { p3 = (pmin < p1) ? p1 : pmin p4 = (pmax < p2) ? pmax : p2 dp = (p4 - p3) / 50 for (p = p3; p <= p4 * 1.001; p += dp) { xyz2(a, b, k1, p, t, p0) print y, z >"'"$TMP1"'" } p1 += PI2 p2 += PI2 print "" >"'"$TMP1"'" } printf "\n\n" >"'"$TMP1"'" } printf "\n\n" >"'"$TMP1"'" # q2 線 for (p = pmin; p <= pmax; p += PI / Np) { t1 = atan2(k1*sqrt(k1^2 - 1)*cos(p), sin(p)) t2 = t1 + PI for (t = t1; t <= t2; t += PI / Ns / 10) { xyz2(a, b, k1, p, t, p0) print y, z >"'"$TMP1"'" } printf "\n\n" >"'"$TMP1"'" } zmin = zmax = 0 # 上端面 for (t = 0; t <= PI2; t += PI / Ns / 2) { xyz2(a, b, k1, pmax, t, p0) print y, z >"'"$TMP1"'" if (z > zmax) zmax = z } printf "\n\n" >"'"$TMP1"'" # 下端面 for (t = 0; t <= PI2; t += PI / Ns / 2) { xyz2(a, b, k1, pmin, t, p0) print y, z >"'"$TMP1"'" if (z < zmin) zmin = z } printf "\n\n" >"'"$TMP1"'" # 輪郭線 for (p = pmin; p <= pmax; p += PI / Np) { t = atan2(k1*sqrt(k1^2 - 1)*cos(p), sin(p)) xyz2(a, b, k1, p, t, p0) print y, z >"'"$TMP1"'" } printf "\n\n" >"'"$TMP1"'" for (p = pmin; p <= pmax; p += PI / Np) { t = atan2(k1*sqrt(k1^2 - 1)*cos(p), sin(p)) xyz2(a, b, k1, p, t + PI, p0) print y, z >"'"$TMP1"'" } # r2 線 # K1 = k1 K1 = -k1 d = (pmax - pmin)/50 for (t0 = 0; t0 <= PI2; t0 += PI / Ns) { pp = 0 for (p = pmin; p <= pmax; p += d) { t = a/b*k1/sqrt(k1^2 - 1)*sqrt(k2^2 - 1)*(p - p0) + t0 xyz2(a, b, K1, p + pc, t, p0) if (p <= pp) { printf "\n" >"'"$TMP2"'" pp = p continue } if (ex <= 0) { if (K1 < 0) print -y, z >"'"$TMP2"'" else print y, z >"'"$TMP2"'" } else printf "\n" >"'"$TMP2"'" pp = p } } exit 0 } END { print "set size ratio 4" print "set nokey" print "set noxtics" print "set noytics" print "set noborder" # print "plot \"'"$TMP1"'\" with lines lt 1, \"'"$TMP2"'\" with lines lt -2" print "plot \"'"$TMP1"'\" with lines lc \"red\", \"'"$TMP2"'\" with lines lc \"black\"" print "" exit 0 } function xyz1(a, k1, p, p0) { x = a * cos(p) y = a * sin(p) #if (k1 < 1) { print "!", k1 >"/dev/stderr"; exit; } z = a / sqrt(k1^2 - 1) * (p - p0) if (k1 < 0) y = -y } function xyz2(a, b, k1, p, t, p0, p1, x1, y1, z1) { p1 = atan2(sqrt(k1^2 - 1), (k1 < 0) ? -1 : 1) xyz1(a, k1, p, p0) x1 = x y1 = y z1 = z ex = cos(p) * cos(t) - cos(p1) * sin(p) * sin(t) ey = sin(p) * cos(t) + cos(p1) * cos(p) * sin(t) ez = -sin(p1) * sin(t) x = x1 - b * ex y = y1 - b * ey z = z1 - b * ez } ' | gnuplot -persist