#author("2024-02-08T14:27:47+09:00","external:moriat","moriat")
#topicpath
** van der Waals の状態方程式について [#y98a85e9]
- コメント~
　気体の状態方程式に限らず、与えられた数式をグラフとして描画できることは、理系として~
ぜひ、身につけておきたい資質である。もちろん、Excel 等を用いて構わない。しかし、世~
の中には大変便利なツールが数多くあるので、それらを利用できると、一層便利である。~
　特に、パラメタに対する依存性を調べたいときには、Excel 等の表計算ソフトより、別の~
ツールの方が便利であることが多い。~
~
　ここでは、ファンデルワールスの状態方程式を例に、式で与えられた情報をグラフにして~
視覚化する方法をあげてみる。~
~
- 目次
-- [[Maximaの場合>#z366c9b4]]
-- [[gnuplot の場合>#w2f9d130]]
-- [[GeoGebra の場合>#d1480be3]]
~
*** ツールの紹介 [#q8e7a89b]
- こんなことができますよ、という紹介
#style(class=table_left)
|Maxima| gnuplot | [[GeoGebra>講義のページ/GeoGebra]] |h
|#youtube_res(4dmvQJBcnIQ)|#youtube_res(kUH134juf_M)|#youtube_res(Yk76i2Trd5A)|
#clear

*** Maxima で計算 [#z366c9b4]
-- 状態方程式;
  
        P(V,RT,N,a,b) := N*RT/( V- b*N) - a * N^2/V^2;
        
:= で関数を定義している。関数の引数は複数でもよい。~
~
--「臨界状態の V, RT, P」 ''を'' 「a, b」 ''で'' 表現;
  
        result1 : solve( [diff(P(V,RT,N,a,b),V)=0, diff(diff(P(V,RT,N,a,b),V),V)=0], [V,RT])[1];
        
solve() の部分で、∂P/∂V=0 と、∂^2P/∂V^2=0 $ \left( \frac{\partial p}{\partial V}=0 , \frac{\partial^2 p}{\partial V^2}=0 \right) $ の両方を満たす V, RT を求めている。~
連立方程式は [ ] 内で、"," でつないで書き連ねる。解く変数も同様に[ ] 内に書く。~
解いた結果は result1 に代入される。: は代入を表す。~
~
  
        P(rhs(result1[1]), rhs(result1[2]), N, a, b);
        
result1 には複数(V, RT)の変数＝値という式が並んでいる。その1番めは result1[1] で取り出し、右辺の値は、rhs() で取り出している。~
それを P() に代入している。~
~
以上の結果を以下にまとめる。
--- 臨界温度 : $ T_c = \frac{8a}{27b R} $
--- 臨界圧力 : $ P_c = \frac{a}{27 b} $
--- その時の体積 : $ V_c = 3 b N $
~
-- 規格化された van der Waals の状態方程式
 
      P( v * 3 * b * N, t * 8 * a / 27 / b / R, N, a, b );
      ratsimp(%);
      factor(%);
 
~
--「臨界状態の V, RT」 ''で'' 「P, a, b」 ''を'' 表現;
  
        result2 : solve( [diff(P(V,RT,N,a,b),V)=0, diff(diff(P(V,RT,N,a,b),V),V)=0], [a,b])[1];
        P(V, RT, N, rhs(result2[1]), rhs(result2[2]));
        

*** gnuplot で描画 [#w2f9d130]
-- gnuplot 用の記述 
        #
        #  van der Waals の状態方程式を
        #  臨界点での温度と体積で規格化する。
        #
        #  パラメタ
        #
        N   = 1.0
        RTc = 1.0
        Vc  = 1.0
        a   = 9.0*RTc*Vc/(8.0*N)
        b   = Vc/(3.0*N)
        set yrange [0:1]
        set xrange [0.4:5]
        
        #
        #  状態方程式 PR:実在気体 PI:理想気体
        #
        PR(V,RT) = N*RT/(V-b*N) - a*N**2/V**2
        PI(V,RT) = N*RT/V
        
        #
        #  確認のため表示
        #  ( gnuplot で整数の割り算には注意。truncation が起こる)
        #
        print 1/3
        print PR(1,1)
        
        #
        #  X11
        #
        set term x11 0
        plot PR(x,2), PR(x,1.5), PR(x,1.2), PR(x,1.1), PR(x,1.0), PR(x,0.90)
        set term x11 1
        plot PI(x,2), PI(x,1.5), PI(x,1.2), PI(x,1.1), PI(x,1.0), PI(x,0.90)
        
        #
        #  gif
        #
        # set term gif small x808080 x909090 xA0A0A0 xB0B0B0 xC0C0C0 xD0D0D0 ; set output "VanDerWaals_PR.gif"
        set term gif small; set output "VanDerWaals_PR.gif"
        plot PR(x,2), PR(x,1.5), PR(x,1.2), PR(x,1.1), PR(x,1.0), PR(x,0.90)
        set term gif small; set output "VanDerWaals_PI.gif"
        plot PI(x,2), PI(x,1.5), PI(x,1.2), PI(x,1.1), PI(x,1.0), PI(x,0.90)
        
        quit
~
-- 結果
理想気体の場合　　　　　　　　　　　　　　　　　　　van der Waals 気体~
&ref(VanDerWaals_PI.gif,,320x240);
&ref(VanDerWaals_PR.gif,,320x240);
~
-- 3次元的な描画
      f(RT,V) = RT/(V-1./3.)-9./8./V/V
      set isosamples 20
      set contour both
      set cntrp le inc 0,0.5,2
      set xrange[0.7:1.5]
      set yrange[0.4:1.5]
      set zrange[0:2]
      splot f(x,y)

*** [[GeoGebra>講義のページ/GeoGebra]] による描画 [#d1480be3]
- Web へ埋め込んだ結果~
#htmlinsert(GeoGebraSample12.txt)
#htmlinsert(GeoGebraSample35.txt)
トップ   編集 差分 添付 複製 名前変更 リロード   新規 検索 最終更新   ヘルプ   最終更新のRSS