問題4.13の解答(回転体の表面積)

の問題4.13の解答。公式をそのまんま利用。

import sympy as sym
from sympy.utilities import lambdify
import scipy as sp
import scipy.integrate

def surface_of_revolution(y, range_):
    dy = sym.diff(y, x)

    y_lambda = lambdify(x, y)
    dy_lambda = lambdify(x, dy)
    def delta_s(x):
        return 2 * sp.pi * y_lambda(x) * sp.sqrt(1 + dy_lambda(x)**2)

    return sp.integrate.quad(delta_s, *range_)

from sympy.abc import x
print surface_of_revolution(x**3 + x, (1., 2.))

(303.5652473801146, 3.370251271174347e-12)

sympyで一般化してあるので、三角錐の表面積とかも求められる。

# a cylinder(S = 2π)
print surface_of_revolution(1, (0., 1.))

# a cone(S = 2π)
print surface_of_revolution(x * sp.sqrt(3) / 3., (0., sp.sqrt(3)))

(6.283185307179586, 6.975736996017264e-14)
(6.283185307179588, 6.975736996017265e-14)