質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう!

新規登録して質問してみよう
ただいま回答率
85.50%
C++

C++はC言語をもとにしてつくられた最もよく使われるマルチパラダイムプログラミング言語の1つです。オブジェクト指向、ジェネリック、命令型など広く対応しており、多目的に使用されています。

Q&A

解決済

1回答

2938閲覧

コンパイルは成功しますが,実行するとAssertion failedといわれます.

chiku_soh

総合スコア12

C++

C++はC言語をもとにしてつくられた最もよく使われるマルチパラダイムプログラミング言語の1つです。オブジェクト指向、ジェネリック、命令型など広く対応しており、多目的に使用されています。

0グッド

0クリップ

投稿2019/01/07 08:44

編集2019/01/07 12:38

Eigen初心者ですが,状態空間方程式を解く際にC++でEigenを用いてプログラムを作成しているのですが,なかなかうまくいきません.解法も含めて教えていただきたいです.
問題は以下の通りですが,

  1. この系の数理模型を作成せよ。

イメージ説明
2. この系の状態空間モデルを作成せよ。ただし、状態x(t) は
以下のように定義せよ。
x(t) = col(y(t), y˙(t))
3. この系のデジタルシミュレータを作成せよ。
4. この系をevent(t0, x0) からevent(t1, x1) に移すための制御入力を見出せ。
ただし、時刻t0 = 0、状態x0 = col(1, 0)、時刻t1 = 1、
状態x1 = col(0, 0) とせよ。

自分での解釈は以下の通りなのですが,
イメージ説明
これが間違ってたら指摘お願いします.

ソースコードは以下の通りです.
なかなかうまくいかなくて困ってます.教えていただけるとありがたいです.

C++

1#include "Eigen/Dense" 2#include "Eigen/Core" 3#include <iostream> 4#include <fstream> 5#include <cmath> 6using namespace Eigen; 7using namespace std; 8 9void RungeKutta(MatrixXd &X, double tt, double dt, MatrixXd A, MatrixXd B, MatrixXd U) { 10 MatrixXd k1 = A*X + B*U; 11 MatrixXd k2 = A*(X + 0.5*k1*dt) + B*U; 12 MatrixXd k3 = A*(X + 0.5*k2*dt) + B*U; 13 MatrixXd k4 = A*(X + k3*dt) + B*U; 14 MatrixXd k = (k1 + 2.0*k2 + 2.0*k3 + k4)*dt / 6.0; 15 X = X + k; 16} 17 18int main() { 19 double m = 1.0; 20 21 MatrixXd A(2, 2); 22 A(0, 0) = 0; 23 A(0, 1) = 0; 24 A(1, 0) = 0; 25 A(1, 1) = 1; 26 MatrixXd B(2, 2); 27 B(0, 0) = 0; 28 B(0, 1) = 0; 29 B(1, 0) = 0; 30 B(1, 1) = 1 / m; 31 32 double dt = 0.01; 33 double tt = 0.0; 34 MatrixXd X(2, 1); 35 X(0, 0) = 1; 36 X(1, 0) = 0; 37 MatrixXd U(2, 1); 38 U(0, 0) = 0; 39 U(1, 0) = -6.25*sin(2.0*tt*M_PI); 40 41 double freq = 1.0; 42 ofstream ofs("outdata.csv"); 43 ofs << "time," << "y" << endl; 44 45 for (int i = 0; i < 1000; i++) { 46 RungeKutta(X, tt, dt, A, B, U); 47 ofs << tt << "," << X(0, 0) << endl; 48 tt += dt; 49 U(0, 1) = -6.25*sin(2.0*tt*M_PI); 50 } 51 52 return 0; 53}

エラー内容は,
Assertion failed: row >= 0 && row < rows() && col >= 0 && col < cols(), file Eigen/src/Core/DenseCoeffsBase.h, line 365
です.
よろしくお願いいたします.

追記!

c++

1#include "Eigen/Dense" 2#include "Eigen/Core" 3#include <iostream> 4#include <fstream> 5#include <cmath> 6using namespace Eigen; 7using namespace std; 8 9 10void RungeKutta(MatrixXd &X, double tt, double dt, MatrixXd A, MatrixXd B, MatrixXd U) { 11 MatrixXd k1 = A*X*dt + -6.25*sin(2.0*tt*M_PI)*B*U*dt; 12 MatrixXd k2 = A*(X + 0.5*k1)*dt + -6.25*sin(2.0*tt*M_PI)*B*U*dt; 13 MatrixXd k3 = A*(X + 0.5*k2)*dt + -6.25*sin(2.0*tt*M_PI)*B*U*dt; 14 MatrixXd k4 = A*(X + k3)*dt + -6.25*sin(2.0*tt*M_PI)*B*U*dt; 15 MatrixXd k = (k1 + 2.0*k2 + 2.0*k3 + k4)*dt / 6.0; 16 X = X + k; 17} 18 19int main() { 20 double m = 1.0; 21 22 MatrixXd A(2, 2); 23 A(0, 0) = 0; 24 A(0, 1) = 1; 25 A(1, 0) = 0; 26 A(1, 1) = 0; 27 MatrixXd B(2, 2); 28 B(0, 0) = 0; 29 B(0, 1) = 0; 30 B(1, 0) = 0; 31 B(1, 1) = 1 / m; 32 33 double dt = 0.001; 34 double tt = 0.0; 35 MatrixXd X(1, 2); 36 X(0, 0) = 1; 37 X(0, 1) = 0; 38 MatrixXd U(1, 2); 39 U(0, 0) = 0; 40 U(0, 1) = 1; 41 42 double freq = 1.0; 43 ofstream ofs("outdata.csv"); 44 ofs << "time," << "x," << "v" << endl; 45 46 for (int i = 0; i <= 1000; i++) { 47 RungeKutta(X, tt, dt, A, B, U); 48 ofs << tt << "," << X(0, 0) << "," << X(0, 1) << endl; 49 tt += dt; 50 } 51 52 return 0; 53}

気になる質問をクリップする

クリップした質問は、後からいつでもMYページで確認できます。

またクリップした質問に回答があった際、通知やメールを受け取ることができます。

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

guest

回答1

0

ベストアンサー

XとUの行と列が逆では?

投稿2019/01/07 12:05

iwanote

総合スコア295

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

chiku_soh

2019/01/07 12:11

すいません。追記し忘れてましたが,逆だと思い試してみましたが,うまくいきませんでした.
iwanote

2019/01/07 12:16

どこでエラーが出るかprintf等で示していただけませんか?
chiku_soh

2019/01/07 12:17

すいません。エラーのほうはStackOverFlowの方で同じような質問があり,解決しました.
chiku_soh

2019/01/07 12:18

実行結果が明らかに軟着陸じゃなく,位置の値がなぜか増加してしまいます.
iwanote

2019/01/07 12:20

では修正後のプログラムと状況を追記してください
iwanote

2019/01/08 08:04

どうおかしくなるか良くわからないので一つづつ数値確かめるしかないかと思います…
chiku_soh

2019/01/09 23:54

ありがとうございました
guest

あなたの回答

tips

太字

斜体

打ち消し線

見出し

引用テキストの挿入

コードの挿入

リンクの挿入

リストの挿入

番号リストの挿入

表の挿入

水平線の挿入

プレビュー

15分調べてもわからないことは
teratailで質問しよう!

ただいまの回答率
85.50%

質問をまとめることで
思考を整理して素早く解決

テンプレート機能で
簡単に質問をまとめる

質問する

関連した質問