物理定数とか

元のエントリ自身とはあんまり関係ないけれども、弾さんのページにておまけの問題があった。

1. 地球の静止軌道の高度を求めなさい。
2. 火星の静止軌道の高度を求めなさい。
3. 火星に軌道エレベーターを作る際に一つ問題となるものがあります。それが何か答えなさい。
4. あなたなら、それをどう解決しますか?

なお、求めた際にはその根拠も述べなさい。数字だけの答えは不合格とします。

簡単そうなのでちょっとやってみる。遠心力=万有引力か。r^3 = GM/(ω^2)。。。即算でr = 42243.141718123246とかやってしまうがこりゃ間違い。これは地球の半径分引いてないからだね。ダメざん orz
そんなことやってたら、google計算機のオフライン版みたいなのがほしくなってくる。オーバーロードとか良く知りもしないくせにpythonでがんばってみる。こんな感じかなー?

"""
base class for constant correction module
"""
from warnings import warn

class C(object):

        def __init__(self, value, unit=""):
                self.value = value
                if type(unit) == str:
                        self.unit = self._str2dict(unit)
                else:
                        self.unit = unit

        def __pos__(self):
                """Inverse only unit'power(for division).
                   this is very ugly . must be fixed"""
                return C(-self.value, 
                                dict([(key, -value) for key, value in self.unit.items()]))

        def __neg__(self):
                return C(-self.value, self.unit)

        def __add__(self, x):
                if not self.unit == x.unit:
                        warn("Constant:: Unit mismatch")
                return C( self.value + x.value, self.unit)

        def __sub__(self, x):
                return self.__add__(-x)

        def __mul__(self, x):
                return C(self.value * x.value, self._merge_unit(self, x))

        def __div__(self, x):
                return C(self.value / x.value, self._merge_unit(self, +x))

        def __str__(self):
                return "%g " % self.value

        def __repr__(self):
                return "%g(%s)" % (self.value, self._dict2str(self.unit))

        @staticmethod
        def _merge_unit(u1, u2):
                """merge two unit-dict"""
                ret = dict(u1.unit)
                for name, pow in u2.unit.items():
                        try:
                                ret[name]+= pow
                        except:
                                ret[name] = pow
                        if ret[name] == 0:
                                del ret[name]
                return ret

        @staticmethod
        def _dict2str(unit_dict):
                """convert dict representation of unit to string"""
                return "".join(["[%s]%d" % pair for pair in sorted(unit_dict.items())])

        @staticmethod
        def _str2dict(unit_str):
                """convert string representation of unit to dict"""
                xs = [(x[0], int(x[1] or 1)) for x in
                                        [x.split(']') for x in 
                                        [x for x in unit_str.split('[')[1:]]]]
                return dict(xs)

キッタない!><

re を使いたくないって理由だけで変態的な分割方法になってしまった。ダメざん。
reで書き直そうかな・・・*1
それと、単位の冪を反転するのに'+'を使ったのはちょっとひどすぎるかも。。後で直そう。

とりあえず、動かしてみよう。

from constant import C

G       = C(6.6742e-11, "[m]3[kg]-1[s]-2")
g       = C(   9.80665, "[m][s]-2")
M_earth = C(5.9742e+24, "[kg]")
M_mars  = C(6.4191e+23, "[kg]")

R_earth = C( 6378.1e+3, "[m]")
R_mars  = C( 3397.0e+3, "[m]")

とか作って

>>> import math
>>> from unit import *
>>> G
6.6742e-11([kg]-1[m]3[s]-2)
>>> g
9.80665([m]1[s]-2)
>>> G*G
4.45449e-21([kg]-2[m]6[s]-4)
>>> G/G
1()

おお。

>>> omega = 2*math.pi /24/60/60
>>> G*M_earth/C(omega**2)
7.53957e+22([m]3[s]-2)

あ、まだ __pow__を実装してなかった・・・

>>> r3 = G*M_earth/C(omega**2)
>>> r = C(r3.value**(1.0/3), "[m]")
>>> r
4.22457e+07([m]1)
>>> (r - R_earth)/C(1000,"[km]-1[m]")
35867.6([km]1)

おおー。ダサいけどできた。
補間機能とあわせて使うと結構快適です。いろいろ作ってみようかな。

        • -

追記:

K       = C( 1/1000.0 , "[k]")
m       = C(   1000.0 , "[milli]")

こんなの追加してみた.
こんな感じ。

>>> omega = C(2*math.pi /24/60/60, "[rad][s]-1")
>>> r3 = G*M_earth/(omega*omega)
>>> r3
7.53957e+22([m]3[rad]-2)
>>> r = C(r3.value**(1.0/3), "[m]")
>>> r
4.22457e+07([m]1)
>>> (r-R_earth)*K
35867.6([k]1[m]1)

意味の無いことにこだわる。てーか [k]1[m]1ってカッコ悪い。 1を消そう。

さらに追記:
きっとこんなの既にあるんだろうなぁ・・・

さらにさらに追記:

>>> omega = C(2*math.pi /24.6/60/60, "[rad][s]-1")
>>> r3 = G*M_mars/(omega*omega)
>>> r = C(r3.value**(1.0/3), "[m]")
>>> (r-R_mars)*K
17020.2([k]1[m]1)

フォボスの高度が

フォボスは、太陽系の衛星の中では最も主星に近い軌道を回っており、その距離は、 火星の表面から 6000 km 以下です。

とのことなので、軌道エレベータが衝突してしまう! ということかな。静止軌道なんだから安心と思ったら公転周期が早いのか。ぶつからないところ無いのかな? 後で調べよう。
#勉強になりました。

*1:こういうところはPerlのほうが好きです