原始人のプログラミング日記

結論から言うと

結果

わー!すごーい!

ソースコード

PhysicalSlime.h

/**
 * @file   PhysicalSlime.h
 * @brief  PhysicalSlime
 * @author Maeda Takumi
 */
#ifndef PHYSICALSLIME_H
#define PHYSICALSLIME_H

#include "main.h"

class PhysicalPoint;

/**
 * @brief  \~english  physical slime
 * @brief  \~japanese 物理的スライム
 * @author \~ Maeda Tkumi
 */
class PhysicalSlime : public KDrawer, private KUpdater {
private:
    /* 中心座標     */ KVector mPosition;
    /* 半径         */ float mRadius;
    /* 水平分割数   */ int mStack;
    /* 垂直分割数   */ int mSlice;
    /* 有効頂点数   */ int mVertexSize;

    /* 質量         */ float mMass;
    /* 頂点毎の質量 */ float mPointMass;

    /* 頂点配列     */ Vector<Vector<KVector>> mVertex;
    /* 法線配列     */ Vector<Vector<KVector>> mNormal;
    /* 頂点毎の質点 */ Vector<Vector<PhysicalPoint*>> mPoints;

    /* 理想体積     */ float mIdealVolume;
    /* 粘度         */ float mViscosity;

    // 体積を計算します。
    float calcVolume() const;
    // 頂点配列のみを移動します(重心の移動に使用します)。
    void translateVertex(const KVector& aPosition);
    // 有効な頂点かどうかを返します。
    bool isValidIndex(const int& aStack, const int& aSlice) const;
    // 有効な垂直分割のインデックスを取得します。
    int getValidSlice(const int& aStack, const int& aSlice) const;
public:
    /**
     * \~english
     * @param aPosition  inital coordinate
     * @param aRadius    inital radius
     * @param aStack     Number of horizonal divisions
     * @param aSlice     Number of vertical divisions
     * @param aMass      mass
     * @param aVisvosity viscosity
     * \~japanese
     * @param aPosition  初期座標
     * @param aRadius    初期半径
     * @param aStack     水平分割数
     * @param aSlice     垂直分割数
     * @param aMass      質量
     * @param aVisvosity 粘度
     */
    PhysicalSlime(
            const KVector& aPosition,
            const float& aRadius,
            const int& aStack,
            const int& aSlice,
            const float& aMass,
            const float& aVisvosity
            );
    virtual ~PhysicalSlime();

    void draw() const override;
    void update() override;

    /**
     * \~english
     * @brief translate center coordinate into argument.
     * @param aPosition new center coordinate
     * \~japanese
     * @brief 中心座標を指定座標に遷移させます。
     * @param aPosition 設定する中心座標
     */
    virtual void translate(const KVector& aPosition);

    /// @brief \~english  move each vertex randomly.
    /// @brief \~japanese 各頂点をランダムに移動します。
    void spike();
    /**
     * \~english
     * @brief apply force.
     * @param aPos   point of application of force
     * @param aForce force vector
     * \~japanese
     * @brief 力を加えます。
     * @param aPos   力の作用点
     * @param aForce 力ベクトル
     */
    void applyForce(
            const KVector& aPos,
            const KVector& aForce
            );
};

#endif /* PHYSICALSLIME_H */

PhysicalSlime.cpp

/**
 * @file   PhysicalSlime.cpp
 * @brief  PhysicalSlime
 * @author Maeda Takumi
 */
#include "PhysicalSlime.h"

#include "PhysicalPoint.h"

PhysicalSlime::PhysicalSlime(
        const KVector& aPosition,
        const float& aRadius,
        const int& aStack,
        const int& aSlice,
        const float& aMass,
        const float& aVisvosity
        ) :
mRadius(aRadius),
mStack(aStack),
mSlice(aSlice),
mVertexSize((mStack - 1) * mSlice + 2),
mMass(aMass),
mPointMass(mMass / mVertexSize),
mVertex(mStack + 1, Vector<KVector>(mSlice)),
mNormal(mStack + 1, Vector<KVector>(mSlice)),
mPoints(mStack + 1, Vector<PhysicalPoint*>(mSlice, NULL)),
mIdealVolume(4.0f / 3.0f * Math::PI * powf(mRadius, 3.0f)),
mViscosity(aVisvosity) {
    const float thetaStack(Math::PI / mStack);
    const float thetaSlice((Math::PI * 2) / mStack);
    // 頂点と法線と質点の割り当て
    for (int i = 0; i <= mStack; ++i) {
        float theta(thetaStack * i);
        float h(mRadius * cos(theta)); // 輪切り高さ
        float r(mRadius * sin(theta)); // 輪切り半径
        for (int j = 0; j < mSlice; ++j) {
            if (!(i == 0 || i == mStack) || j == 0) { // 重なる頂点は作らない
                float phi(thetaSlice * j);
                KVector vertex(r * sin(phi), h, r * cos(phi));
                mVertex[i][j] = vertex;
                mNormal[i][j] = vertex.normalization();
                mPoints[i][j] = new PhysicalPoint(mPointMass, vertex);
            }
        }
    }
    translate(aPosition);

    spike();
}

PhysicalSlime::~PhysicalSlime() {
    for (Vector<PhysicalPoint*>& i : mPoints) {
        for (PhysicalPoint*& j : i) {
            if (j) delete j;
        }
    }
}

void PhysicalSlime::draw() const {
    glColor4ub(0xfc, 0xe2, 0xc4, 0x11);

    glBegin(GL_TRIANGLE_STRIP);
    for (int i = 0; i < mStack; ++i) {
        for (int j = 0; j <= mSlice; ++j) {
            int slice(getValidSlice(i, j));
            KVector n(mNormal[i][slice]);
            KVector v(mPoints[i][slice]->position());
            glNormal3f(DEPLOY_VEC(n));
            glVertex3f(DEPLOY_VEC(v));

            slice = getValidSlice(i + 1, j);
            n = mNormal[i + 1][slice];
            v = mPoints[i + 1][slice]->position();
            glNormal3f(DEPLOY_VEC(n));
            glVertex3f(DEPLOY_VEC(v));
        }
    }
    glEnd();
}

void PhysicalSlime::update() {
    const float cVol(calcVolume());

    float correction((mIdealVolume - cVol) / mIdealVolume);
    correction = Math::min(Math::max(-1.0f, correction), 1.0f) / mVertexSize;
    println(correction * mVertexSize);

    KVector centroid;
    for (int i = 0; i <= mStack; ++i) {
        for (int j = 0; j < mSlice; ++j) {
            if (isValidIndex(i, j)) {
                KVector position(mPoints[i][j]->position());

                mNormal[i][j] = (position - mPosition).normalization();
                mPoints[i][j]->applyForce(mNormal[i][j] * mRadius * correction);

                KVector resilience(mVertex[i][j] - position);
                mPoints[i][j]->applyForce(resilience / sqrtf(mVertexSize));

                // 速度を0にする力(柔らかさを変更できるようにする)
                mPoints[i][j]->applyForce(-mPoints[i][j]->velocity() * mPointMass * 1.0_s * mViscosity);

                centroid += position;
            }
        }
    }

    // 重心の更新
    translateVertex(centroid / mVertexSize);
}

float PhysicalSlime::calcVolume() const {
    const KVector y_axis(mVertex[mStack][0] - mVertex[0][0]);

    float volume(0.0f);

    // 上面・底面(輪切りにして計算(積分のイメージ))
    float tArea(0.0f), bArea(0.0f);
    KVector tCenter = mVertex[mStack][0], bCenter;
    for (int i = 1; i <= mStack; ++i) {
        bArea = 0.0f;
        bCenter = KVector();
        if (i != mStack) {
            for (int j = 0; j < mSlice; ++j) bCenter += mPoints[i][getValidSlice(i, j)]->position();
            bCenter /= mSlice;

            KVector point, next = mPoints[i][0]->position();
            for (int j = 0; j < mSlice; ++j) {
                point = next;
                next = mPoints[i][getValidSlice(i, j)]->position();
                KVector a(point - bCenter), b(next - bCenter); // 中心から各頂点までのベクトル
                bArea += a.length() * b.length() * a.angle(b);
            }
        }
        // (半径^2 * θ) / 2(扇形の面積)
        bArea /= 2; // まとめて割る

        KVector heightV((bCenter - tCenter).extractParallel(y_axis));
        float aveArea((tArea + bArea) / (tArea && bArea ? 2 : 3));
        volume += Math::sign(heightV.dot(y_axis)) * aveArea * heightV.length();

        tArea = bArea;
        tCenter = bCenter;
    }
    return volume;
}

void PhysicalSlime::translateVertex(const KVector& aPosition) {
    const KVector move(aPosition - mPosition);
    mPosition += move;
    for (Vector<KVector>& i : mVertex) {
        for (KVector& j : i) {
            j += move;
        }
    }
}

bool PhysicalSlime::isValidIndex(const int& aStack, const int& aSlice) const {
    return mPoints[aStack][aSlice];
}

int PhysicalSlime::getValidSlice(const int& aStack, const int& aSlice) const {
    if (aSlice == mSlice || !isValidIndex(aStack, aSlice)) return 0;
    return aSlice;
}

void PhysicalSlime::translate(const KVector & aPosition) {
    const KVector move(aPosition - mPosition);
    mPosition += move;
    for (int i = 0; i <= mStack; ++i) {
        for (int j = 0; j < mSlice; ++j) {
            if (isValidIndex(i, j)) {
                mVertex[i][j] += move;
                mPoints[i][j]->translate(mPoints[i][j]->position() + move);
            }
        }
    }
}

void PhysicalSlime::spike() {
    for (int i = 0; i <= mStack; ++i) {
        for (int j = 0; j < mSlice; ++j) {
            if (isValidIndex(i, j)) {
                KVector randVec(random(300) - 150, random(300) - 150, random(300) - 150);
                mPoints[i][j]->translate(mPoints[i][j]->position() + randVec.normalization() * mRadius * 0.5);
            }
        }
    }
}

void PhysicalSlime::applyForce(const KVector& aPos, const KVector& aForce) {
    // 力の作用店に最も近い頂点を探す。
    PhysicalPoint * nearPoint(NULL);
    float nearLen(0xffffff);
    for (int i = 0; i <= mStack; ++i) {
        for (int j = 0; j < mSlice; ++j) {
            if (isValidIndex(i, j)) {
                float len((mPoints[i][j]->position() - aPos).length());
                if (len < nearLen) {
                    nearLen = len;
                    nearPoint = mPoints[i][j];
                }
            }
        }
    }
    if (nearPoint) { // 力の届く範囲内
        // 最近頂点からの距離により伝わる力を調整する。
        // 粘度が強いほど、力を伝わりやすくする。
        KVector nearPos(nearPoint->position());
        float diameter(mRadius * 2.0f);
        for (int i = 0; i <= mStack; ++i) {
            for (int j = 0; j < mSlice; ++j) {
                if (isValidIndex(i, j)) {
                    float len((mPoints[i][j]->position() - nearPos).length());
                    float correction((Math::max((diameter - len) / diameter, 0.0f) + mViscosity) / mVertexSize);
                    mPoints[i][j]->applyForce(-(nearPos - aPos).normalization() * aForce.length() * correction);
                }
            }
        }
    }
}

解説

注意

このプログラムはしっかりとした物理法則にしたがうものではなく、学術的な理論とは程遠い勘によって、それなりにスライムらしい挙動を実装しただけです。

概要

このプログラムの要点は、毎フレーム各頂点に3つの力を加えていることです(実質それだけです)。

体積を一定にしようとする力
球体になろうとする力
慣性を打ち消す力

1つめの力は、理想体積と計算で求める体積との比率により大きさが変動する力です。
2つめの力は、現在の頂点位置から理想的な球形への頂点位置と戻ろうとする力です。
3つめの力は、球形を安定させる力です。

1と2の力は物理的に理解できるのですが、3の力が物理法則に従っているかどうかが分かりません。
かといって3の力がない場合、球形が安定せず体積が無限に大きく(小さく)なってしまいます。

終わりに

試行錯誤が多く完成したときはテンションが異常に上がったのですが、いざ解説を書こうとすると大変でテンション下がりました。
解説はおざなりでしたが、パラメータや加える力の式を変えることで大きく挙動が変わると思うので今後とも作りこんでいこうと思います。

嘘です。最弱です。

文明が成熟した現代において、素のC++で書かれたライブラリは作者のひいき目から見ても弱いです。
C++自体は弱くないので好きです。

Kawazuiniライブラリ

GitHubでソースコードを公開しています。

github.com

特徴

使うにはKApplicationを継承したクラスのupdate()とdraw()を実装するだけ。
ライブラリのクラスはすべて「K-」ではじまる。
WindowsAPIが使える。
OpenGLが使える。
GLUTを使わずに作成されている。

使用方法

前提として

WindowsAPIを使用しているのでWindows限定です。
C++11に準拠して書いているつもりです*1。
ライブラリを使用して起きるいかなる不具合について責任は取りません。

事前準備

リンクの必要なライブラリ

lgdi32
lgdiplus
lOle32
lmsimg32
lwinmm
lOpenGL32
lGLU32

遠い昔の記憶によれば、MinGWのC++環境にすべて入っているはず。
なので、特別にライブラリをインストールする必要はないと思います。

具体的なソースコード

main.cpp

#include "Kawazuini.h"
#include "MyApp.h"

int WINAPI _tWinMain(
        HINSTANCE aInst,
        HINSTANCE aPrevInst,
        LPSTR aCmdLine,
        int aCmdShow
        ) {
    KWindow::MainArgs args{aInst, aPrevInst, aCmdLine, aCmdShow};

    try {
        KWindow window(args, KRect(800, 800), "ウィンドウタイトル", true);
        MyApp app(window);
        window.show();
        app.start(50);
    } catch (Error& e) {
        message(e.what(), "初期化エラー");
        return 1;
    };

    return 0;
}

MyApp.h

#ifndef MYAPP_H 
#define MYAPP_H 

#include "Kawazuini.h"

class MyApp final : public KApplication {
private:
    KGLTextLabel mFPSLabel;
public:
    MyApp(KWindow& aWindow);
    ~MyApp() = default;

    void draw() const override;
    void update() override;
};

#endif /* MYAPP_H */

MyApp.cpp

static const KImage IMG_CHARSET_MEDIUM(IDI_CHARSET_MEDIUM);
static const KCharset CHARSET_MEDIUM(IMG_CHARSET_MEDIUM, 16);

MyApp::MyApp(KWindow& aWindow) :
KApplication(aWindow),
mFPSLabel(KVector(), CHARSET_MEDIUM, loadString(IDS_FPS) + "100") {
    mFrontUI.addContent(mFPSLabel);
}

void MyApp::draw() const {
    KDrawer::DRAW();
}

void MyApp::update() {
    static const String TEXT_FPS("FPS : ");
    mFPSLabel.setText(TEXT_FPS + toString(mFrame));
}

これだけでウィンドウが生成されて画面上にフレームレートが描画されるはずです。*2 *3

終わりに

ライブラリは趣味で作っているだけなので、温かい目で見ていただければいいと思います。
また今回説明していないことも、GitHubにあげているドキュメントに書いてあったりするので、もし使おうと思ったらそちらを見てください。
適当でしたが、こんな感じでやっていきます。

*1: MinGWでのコンパイル以外は検証していません。

*2:OpenGLでUIが作成されているため、文字列描画のために文字列描画用の画像が必要になります。

*3:ソースコードはコンパイルしてないので失敗するかもしれません。

原始人のプログラミング日記

println(++C++); // println(C++), ++C;

C++でスライムを作ってみた

結論から言うと

結果

ソースコード

解説

注意

概要

終わりに

ぼくのかんがえたさいきょうの　らいぶらり

嘘です。最弱です。

Kawazuiniライブラリ

GitHubでソースコードを公開しています。

特徴

使用方法

前提として

事前準備

具体的なソースコード

終わりに