MATLAB,Python,Scilab,Julia比較 その55【状態空間モデル⑬】

バックナンバーはこちら。
https://www.simulationroom999.com/blog/compare-matlabpythonscilabjulia-backnumber/

はじめに
登場人物
【再掲】ブロック図＆漸化式
Cコード化
シミュレーション結果
まとめ

はじめに

前回から、運動方程式を元に、

ブロック図を作成
上記ブロック図の離散化
上記を元に漸化式導出

漸化式まで出てるとCコード化が可能になる。
このCコードを今回確認する。

登場人物

博識フクロウのフクさん

イラストACにて公開の「kino_k」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=iKciwKA9&area=1

エンジニア歴8年の太郎くん

イラストACにて公開の「しのみ」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=uCKphAW2&area=1

【再掲】ブロック図＆漸化式

太郎くん

今回は運動方程式のCコードを見せてくれるんだっけ？

フクさん

そうだね。
まずは各ブロック図、漸化式を再掲しておこう。

運動方程式のブロック図

運動方程式のブロック図(離散化)

速度の漸化式
\(\boldsymbol{\dot{v}}(t)=\displaystyle \frac{1}{m}\boldsymbol{F}(t)\)
\(\boldsymbol{v}(t)=\boldsymbol{\dot{v}}(t)\Delta T + \boldsymbol{v}(t-\Delta T) \)

距離の漸化式
\(\boldsymbol{\dot{s}}(t)=\boldsymbol{v}(t)\)
\(\boldsymbol{s}(t)=\boldsymbol{\dot{s}}(t)\Delta T + \boldsymbol{s}(t-\Delta T) \)

Cコード化

フクさん

上記の情報を元に作成したCコードが以下になる。

#include <stdio.h>

double dotV;
double V;
double Vzi;
double dotS;
double S;
double Szi;

const double m=1;

int ObjectMovementModel(double F, double deltaT)
{
	int ret = -1;
	if( m != 0 ){
		dotV = 1/m * F;
		V = dotV * deltaT + Vzi;
		Vzi = V;
		
		dotS = V;
		
		S = dotS * deltaT + Szi;
		Szi = S;
	}
	
	return ret;
}

int main()
{
	int i;
	double F = 0;
	printf("F,v,s\n");
	for( i = 0; i < 10000; i++ ){
		if( i >= 5000 ){
			F = 1;
		}
		ObjectMovementModel( F, 0.001 );
		printf("%lf,%lf,%lf\n",F, V, S);
	}
	
	return 0;
}