LoginSignup
20
10

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 1 year has passed since last update.

@tommyecguitar(とみー)

マクスウェルの悪魔ゲーム

Posted at

まずマクスウェルの悪魔とは

マクスウェルの悪魔とは、あのマクスウェルが提唱した、熱力学第二法則に関する思考実験のタイトルもしくはそこに登場する架空の存在です。

  1. 均一な温度の気体で満たされた容器を用意する。 このとき温度は均一でも個々の分子の速度は決して均一ではないことに注意する。
  2. この容器を小さな穴の空いた仕切りで2つの部分 A, B に分離し、個々の分子を見ることのできる「存在」がいて、この穴を開け閉めできるとする。
  3. この存在は、素早い分子のみを A から B へ、遅い分子のみを B から A へ通り抜けさせるように、この穴を開閉するのだとする。
  4. この過程を繰り返すことにより、この存在は力学的意味で仕事をすることなしに、 A の温度を下げ、 B の温度を上げることができる。 これは熱力学第二法則と矛盾する。

悪魔になって部屋のエントロピーを下げるゲームができそうだなと思って作ってみました。
こんなこと他の誰かも思いついてるんじゃないかと思いましたが、ふと浮かんだのでやってみます。

ルール

  1. 舞台は真ん中を壁で仕切られた、断熱・定圧の箱の中です。
  2. あなたは悪魔です。
  3. あなたは、箱の真ん中の壁の真ん中にいます。
  4. 壁の真ん中には小さな扉があります。
  5. あなたはその扉を開け閉めできます。
  6. 箱の中は高温で運動の激しい粒子と、低温で運動の穏やかな粒子が動き回っています。
  7. 扉を開けたり閉めたりして、高温の粒子を右の部屋に移動させ、低温の粒子を左の部屋に分離したらクリアです。
  8. 最後は扉を閉じてください。

MATLABでの設計

コードは公開しています。

基本的にハンドルクラスを継承して、クラスを関数に入力したときに参照渡しされるようにします。

箱

箱のパラメータや状態を格納します。
箱のメソッドは非常に簡単で、パラメータの設定と、ドアの開け閉めフラグの反転のみです。

コード
fieldBox.m
classdef fieldBox < handle
    % Box クラス
    % 箱の状態とパラメータを格納します。
    % 箱と粒子の衝突は粒子クラスでやります。

    properties(SetAccess=private)
        width (1,1) single = 100 % 箱の横サイズ
        height (1,1) single = 100 % 箱の縦サイズ
        wallXPos (1,1) single = 50 % 壁の横位置
        doorYPos (1,1) single = 50 % 扉の縦位置
        doorHeight (1,1) single = 10 % 扉の高さ
        isOpenDoor (1,1) logical = false % 扉が開いているかどうか
    end

    methods(Access=public)
        function obj = fieldBox(width,height,wallXPos,doorYPos,doorHeight)
            % このクラスのインスタンスを作成
            % 現在は箱サイズと扉高さのみの指定
            obj.width = width;
            obj.height = height;
            obj.doorHeight = doorHeight;

            obj.wallXPos = wallXPos;
            obj.doorYPos = doorYPos;
            obj.isOpenDoor = false;

        end

        function DoorSwitch(obj)
            obj.isOpenDoor = ~obj.isOpenDoor;
        end
    end
end

粒子

粒子は数が多いので絶対にハンドルを継承した方がいいです。

工夫点としては、行列計算の得意なMATLABの特性を活かすため、Array of structureではなくStructure of arrayの形式にして、反射の計算などをいっぺんにできるようにしています。

また、単純に入射角=反射角の形式で反射を計算すると、いつまでも扉の方に飛んでこない粒子ができるので、反射した直後に速度方向をずらす工夫をしています。

コード
particles.m
classdef particles < handle
    % 粒子クラス
    % 粒子の状態とパラメータを保管する。

    properties(SetAccess=private)
        velocity (:,2) single
        position (:,2) single
        loParticleNum (1,1) int32
        hiParticleNum (1,1) int32
        speed (:,1) single
    end

    methods
        function obj = particles(hiParticleNum, loParticleNum, speedLo, speedHi, iBox)
            arguments
                hiParticleNum (1,1) int32
                loParticleNum (1,1) int32
                speedLo (1,1) single
                speedHi (1,1) single
                iBox (1,1) fieldBox
            end
            
            obj.loParticleNum = loParticleNum;
            obj.hiParticleNum = hiParticleNum;
            totalNum = hiParticleNum + loParticleNum;

            % 粒子を生成する。
            angLo = rand(obj.loParticleNum,1) * 2 * pi;
            angHi = rand(obj.hiParticleNum,1) * 2 * pi;
            velLo = [cos(angLo), sin(angLo)] .* speedLo;
            velHi = [cos(angHi), sin(angHi)] .* speedHi;
            obj.velocity = [velLo; velHi];
            obj.speed = [zeros(obj.loParticleNum,1) + speedLo;
                zeros(obj.hiParticleNum,1) + speedHi];
            
            % 箱の領域内の点の位置を計算する.箱の左下が原点。
            obj.position = (rand(totalNum,2)*0.99 + 0.005) .* [iBox.width, iBox.height];
        end

        function bound(obj,iBox,dt)
            arguments
                obj (1,1) particles
                iBox (1,1) fieldBox
                dt (1,1) single
            end
            
            % 速度に応じて移動
            diffPos = obj.velocity * dt;
            obj.position = obj.position + diffPos;

            % 右の壁から外に行ったやつを反射
            isRightBound = obj.position(:,1) >= iBox.width;
            obj.position(isRightBound,1) = -obj.position(isRightBound,1) + 2 * iBox.width;

            % 左の壁から外に行ったやつを反射
            isLeftBound = obj.position(:,1) <= 0;
            obj.position(isLeftBound,1) = -obj.position(isLeftBound,1);

            % 上の壁から外に行ったやつを反射
            isRoofBound = obj.position(:,2) >= iBox.height;
            obj.position(isRoofBound,2) = -obj.position(isRoofBound,2) + 2 * iBox.height;

            % 下の壁から外に行ったやつを反射
            isFloorBound = obj.position(:,2) <= 0;
            obj.position(isFloorBound,2) = -obj.position(isFloorBound,2);

            % 真ん中の壁にヒットしたやつ
            xFromWall = obj.position(:,1) - iBox.wallXPos;
            isWallBound = xFromWall .* (xFromWall - diffPos(:,1)) <= 0;
            
            % ドアが開いてたらドア通り抜けるので反射判定をオフる。
            if iBox.isOpenDoor
                % ドアを通り抜けないやつを計算
                posOnTheWall = obj.position(:,2) - diffPos(:,2) .*...
                    ((obj.position(:,1) - iBox.wallXPos) ./ diffPos(:,1));
                isNotThroughDoor = abs(posOnTheWall-iBox.doorYPos) > iBox.doorHeight/2;
                isWallBound = isWallBound & isNotThroughDoor;
            end

            % 壁ヒット反射
            obj.position(isWallBound,1) = -obj.position(isWallBound,1) + 2 * iBox.wallXPos;

            % 水平方向の反射の速度更新
            % 反射方向のせいで永遠にドアに来ない粒子が出ないように、進行方向をずらす。
            isHRef = isRightBound | isLeftBound | isWallBound;
            newVel = obj.velocity(isHRef,:) .* [-1,1];
            newVel = newVel + obj.speed(isHRef) .* (2*single(newVel >= 0)-1) .* rand(size(newVel,1),2);
            newVel = obj.speed(isHRef) .* newVel ./ vecnorm(newVel,2,2);
            obj.velocity(isHRef,:) = newVel;

            % 垂直方向の反射の速度更新
            isVRef = isRoofBound | isFloorBound;
            newVel = obj.velocity(isVRef,:) .* [1,-1];
            newVel = newVel + obj.speed(isVRef) .* (2*single(newVel >= 0)-1) .* rand(size(newVel,1),2);
            newVel = obj.speed(isVRef) .* newVel ./ vecnorm(newVel,2,2);
            obj.velocity(isVRef,:) = newVel;
        end
    end
end

ゲームクラス

ゲームの状態を格納して、クリアの判定などをします。

ゲームの初期設定と、実行関数があります。
ウィンドウコールバックを使って扉の開け閉めを設計しています。

コード
game.m
classdef MaxwellDemonGame < handle
    %UNTITLED このクラスの概要をここに記述
    %   詳細説明をここに記述

    properties
        lBox
        lParticle
        fig
        doorPlot
        pScatter
        figStatus
        gameStatus
    end

    methods(Access=private)
        function SetFigure(obj)
            obj.fig = figure;
            obj.fig.WindowKeyPressFcn = @obj.figureCallBack;
            obj.fig.CloseRequestFcn = @obj.figureClosed;
            ax = gca(obj.fig);
            rectangle(ax,"Position",[0,0,obj.lBox.width,obj.lBox.height],"LineWidth",2,"EdgeColor","k");
            xlim(ax,[0,obj.lBox.width])
            ylim(ax,[0,obj.lBox.height])
            axis(ax,"equal")
            hold(ax,'on')
            plot(ax,[obj.lBox.wallXPos,obj.lBox.wallXPos],[obj.lBox.height,obj.lBox.doorYPos+obj.lBox.doorHeight/2],...
                "LineWidth",2,"Color","k")
            plot(ax,[obj.lBox.wallXPos,obj.lBox.wallXPos],[0,obj.lBox.doorYPos-obj.lBox.doorHeight/2],...
                "LineWidth",2,"Color","k")
            obj.doorPlot = plot(ax,[obj.lBox.wallXPos,obj.lBox.wallXPos],[obj.lBox.doorYPos-obj.lBox.doorHeight/2,obj.lBox.doorYPos+obj.lBox.doorHeight/2],...
                "LineWidth",2,"Color","r");
            obj.pScatter = scatter(ax,obj.lParticle.position(:,1), obj.lParticle.position(:,2), 20, vecnorm(obj.lParticle.velocity,2,2), "filled");
            colormap(ax,"flag")
            axis off
        end

                function judge(obj)
            loPosJudge = obj.lParticle.position(1:obj.lParticle.loParticleNum,1) < obj.lBox.wallXPos;
            hiPosJudge = obj.lParticle.position(obj.lParticle.loParticleNum+1:obj.lParticle.loParticleNum+obj.lParticle.hiParticleNum,1) < obj.lBox.wallXPos;
            loJudge = nnz(hiPosJudge) == 0 & nnz(loPosJudge) == obj.lParticle.loParticleNum;
            hiJudge = nnz(loPosJudge) == 0 & nnz(hiPosJudge) == obj.lParticle.hiParticleNum;
            obj.gameStatus = and(or(loJudge,hiJudge),~obj.lBox.isOpenDoor);
        end

        function figureCallBack(obj,~,event)
            if event.Character == "s"
                obj.lBox.DoorSwitch();
                obj.doorPlot.Visible = ~obj.doorPlot.Visible;
            end
        end

        function figureClosed(obj,~,~)
            obj.figStatus = false;
            delete(gcf)
        end
    end

    methods
        function obj = MaxwellDemonGame(iBox, iParticle)
            obj.lBox = iBox;
            obj.lParticle = iParticle;
            
            obj.SetFigure();

            obj.figStatus = true;
            obj.gameStatus = false;
        end

        function res = run(obj)
            
            dt = 1/60;

            while obj.figStatus && ~obj.gameStatus
                preTime = cputime;
                obj.lParticle.bound(obj.lBox, dt);
                obj.pScatter.XData = obj.lParticle.position(:,1);
                obj.pScatter.YData = obj.lParticle.position(:,2);
                drawnow
                % 判定
                obj.judge();
                dt = cputime - preTime;
            end
            
            if obj.gameStatus
                ax = gca(obj.fig);
                text(ax,obj.lBox.width/2,obj.lBox.height/2,"CLEAR!! You Are Demon!!","FontSize",30,"Color","r","HorizontalAlignment","center","VerticalAlignment","middle");
            end

            res = true;
        end
    end
end

まとめ

ネタ記事&ネタコードなので冷やかし程度に遊んでみてください。
粒子の数や速度をいじって激ムズゲームにすることもできるので、各自工夫して遊んでください。

20
10
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
20
10

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?