@@ -48,7 +48,13 @@ void evaluateResidual(std::vector<double> y_, std::vector<double> yp_, std::vect
4848
4949 constexpr bool ref_frame_ = true ;
5050
51- double omega = wb_ - mp_ * y_[4 ];
51+ using ScalarT = double ;
52+
53+ ScalarT omega = wb_ - mp_ * y_[3 ];
54+
55+ ScalarT vbd_in = std::cos (y_[15 ]) * y_[1 ] + std::sin (y_[15 ]) * y_[2 ];
56+ ScalarT vbq_in = -std::sin (y_[15 ]) * y_[1 ] + std::cos (y_[15 ]) * y_[2 ];
57+
5258 if (ref_frame_)
5359 {
5460 f_[0 ] = omega - y_[0 ];
@@ -58,39 +64,47 @@ void evaluateResidual(std::vector<double> y_, std::vector<double> yp_, std::vect
5864 f_[0 ] = 0.0 ;
5965 }
6066
61- f_[1 ] = cos (y_[3 ]) * y_[14 ] - sin (y_[3 ]) * y_[15 ];
62- f_[2 ] = sin (y_[3 ]) * y_[14 ] + cos (y_[3 ]) * y_[15 ];
67+ // output
68+ // current transformed to common frame
69+ f_[1 ] = std::cos (y_[15 ]) * y_[13 ] - std::sin (y_[15 ]) * y_[14 ];
70+ f_[2 ] = std::sin (y_[15 ]) * y_[13 ] + std::cos (y_[15 ]) * y_[14 ];
6371
64- double vbd_in = cos (y_[3 ]) * y_[1 ] + sin (y_[3 ]) * y_[2 ];
65- double vbq_in = -sin (y_[3 ]) * y_[1 ] + cos (y_[3 ]) * y_[2 ];
72+ // ### Internal Componenets ##
73+ // P and Q equations
74+ f_[3 ] = -yp_[3 ] + wc_ * (y_[11 ] * y_[13 ] + y_[12 ] * y_[14 ] - y_[3 ]);
75+ f_[4 ] = -yp_[4 ] + wc_ * (-y_[11 ] * y_[14 ] + y_[12 ] * y_[13 ] - y_[4 ]);
6676
67- f_[ 3 ] = -yp_[ 3 ] + omega - y_[ 0 ];
68- f_[ 4 ] = -yp_[ 4 ] + wc_ * ( y_[12 ] * y_[ 14 ] + y_[ 13 ] * y_[ 15 ] - y_[ 4 ]) ;
69- f_[ 5 ] = -yp_[ 5 ] + wc_ * (-y_[ 12 ] * y_[ 15 ] + y_[ 13 ] * y_[ 14 ] - y_[ 5 ] );
77+ // Voltage control
78+ ScalarT vod_star = Vn_ - nq_ * y_[4 ] ;
79+ ScalarT voq_star = static_cast <ScalarT>( 0.0 );
7080
71- double vod_star = Vn_ - nq_ * y_[5 ];
72- double voq_star = 0.0 ;
81+ f_[ 5 ] = -yp_[ 5 ] + vod_star - y_[11 ];
82+ f_[ 6 ] = -yp_[ 6 ] + voq_star - y_[ 12 ] ;
7383
74- f_[ 6 ] = -yp_[ 6 ] + vod_star - y_[12 ];
75- f_[ 7 ] = -yp_[ 7 ] + voq_star - y_[13 ];
84+ ScalarT ild_star = F_ * y_[ 13 ] - wb_ * Cf_ * y_[ 12 ] + Kpv_ * ( vod_star - y_[11 ]) + Kiv_ * y_[ 5 ];
85+ ScalarT ilq_star = F_ * y_[ 14 ] + wb_ * Cf_ * y_[ 11 ] + Kpv_ * ( voq_star - y_[12 ]) + Kiv_ * y_[ 6 ];
7686
77- double ild_star = F_ * y_[14 ] - wb_ * Cf_ * y_[13 ] + Kpv_ * (vod_star - y_[12 ]) + Kiv_ * y_[6 ];
78- double ilq_star = F_ * y_[15 ] + wb_ * Cf_ * y_[12 ] + Kpv_ * (voq_star - y_[13 ]) + Kiv_ * y_[7 ];
87+ // Current control
88+ f_[7 ] = -yp_[7 ] + ild_star - y_[9 ];
89+ f_[8 ] = -yp_[8 ] + ilq_star - y_[10 ];
7990
80- f_[ 8 ] = -yp_[ 8 ] + ild_star - y_[10 ];
81- f_[ 9 ] = -yp_ [9 ] + ilq_star - y_[11 ];
91+ ScalarT vid_star = -wb_ * Lf_ * y_[ 10 ] + Kpc_ * ( ild_star - y_[9 ]) + Kic_ * y_[ 7 ];
92+ ScalarT viq_star = wb_ * Lf_ * y_ [9 ] + Kpc_ * ( ilq_star - y_[10 ]) + Kic_ * y_[ 8 ];
8293
83- double vid_star = -wb_ * Lf_ * y_[11 ] + Kpc_ * (ild_star - y_[10 ]) + Kic_ * y_[8 ];
84- double viq_star = wb_ * Lf_ * y_[10 ] + Kpc_ * (ilq_star - y_[11 ]) + Kic_ * y_[9 ];
94+ // Output LC Filter
95+ f_[9 ] = -yp_[9 ] - (rLf_ / Lf_) * y_[9 ] + omega * y_[10 ] + (vid_star - y_[11 ]) / Lf_;
96+ f_[10 ] = -yp_[10 ] - (rLf_ / Lf_) * y_[10 ] - omega * y_[9 ] + (viq_star - y_[12 ]) / Lf_;
8597
86- f_[10 ] = -yp_[10 ] - (rLf_ / Lf_) * y_[10 ] + omega * y_[11 ] + (vid_star - y_[12 ]) / Lf_ ;
87- f_[11 ] = -yp_[11 ] - (rLf_ / Lf_) * y_[11 ] - omega * y_[10 ] + (viq_star - y_[13 ]) / Lf_ ;
98+ f_[11 ] = -yp_[11 ] + omega * y_[12 ] + ( y_[9 ] - y_[13 ]) / Cf_ ;
99+ f_[12 ] = -yp_[12 ] - omega * y_[11 ] + ( y_[10 ] - y_[14 ]) / Cf_ ;
88100
89- f_[12 ] = -yp_[12 ] + omega * y_[13 ] + (y_[10 ] - y_[14 ]) / Cf_;
90- f_[13 ] = -yp_[13 ] - omega * y_[12 ] + (y_[11 ] - y_[15 ]) / Cf_;
101+ // Output Connector
102+ f_[13 ] = -yp_[13 ] - (rLc_ / Lc_) * y_[13 ] + omega * y_[14 ] + (y_[11 ] - vbd_in) / Lc_;
103+ f_[14 ] = -yp_[14 ] - (rLc_ / Lc_) * y_[14 ] - omega * y_[13 ] + (y_[12 ] - vbq_in) / Lc_;
91104
92- f_[14 ] = -yp_[14 ] - (rLc_ / Lc_) * y_[14 ] + omega * y_[15 ] + (y_[12 ] - vbd_in) / Lc_;
93- f_[15 ] = -yp_[15 ] - (rLc_ / Lc_) * y_[15 ] - omega * y_[14 ] + (y_[13 ] - vbq_in) / Lc_;
105+ // Rotor difference angle
106+ if (!ref_frame_)
107+ f_[15 ] = -yp_[15 ] + omega - y_[0 ];
94108}
95109
96110// Function that computes the Jacobian via automatic differentiation
0 commit comments